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Abstract 



Hysteresis loops are often seen in experiments at first order phase transforma- 



■ tions, when the system goes out of equilibrium. They may have a macroscopic 

o 

jump (roughly as in the supercooling of liquids) or they may be smoothly 
varying (as seen in most magnets). We have studied the nonequilibrium zero- 

a 

i . temperature random-field Ising-model as a model for hysteretic behavior at 

. 

first order phase transformations. As disorder is added, one finds a transition 

O 

where the jump in the magnetization (corresponding to an infinite avalanche) 
^ . decreases to zero. At this transition we find a diverging length scale, power 

law distributions of noise (avalanches) and universal behavior. We expand the 
critical exponents about mean- field theory in 6 — e dimensions. Using a map- 
ping to the pure Ising model, we Borel sum the 6 — e expansion to 0(e 5 ) for the 
correlation length exponent. We have developed a new method for directly 
calculating avalanche distribution exponents, which we perform to O(e). Nu- 
merical exponents in three, four, and five dimensions are in good agreement 
with the analytical predictions. Some suggestions for further analyses and 
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experiments are also discussed. 
PACS numbers: 75.60.Ej, 64.60.Ak, 81.30.Kf 
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I. INTRODUCTION 



The modern field of disordered systems has its roots in dirt. An important effect of 
disorder is the slow relaxation to equilibrium seen in many experimental systems [1]. This 
paper is an attempt to unearth universal, nonequilibrium collective behavior buried in the 
muddy details of real materials and inherently due to their tendency to remain far from 
equilibrium on experimental time scales. In particular, we focus on two distinctly nonequi- 
librium effects: (a) the avalanche response to an external driving force and (b) the internal 
history dependence of the system (hysteresis). 

Systems far from equilibrium often show interesting memory effects not present in equi- 
librium systems. Far from equilibrium, the system will usually occupy some metastable 
state that has been selected according to the history of the system. Jumps over large free 
energy barriers to reach a more favorable state are unlikely. The system will move through 
the most easily accessible local minima in the free energy landscape as an external driving 
field is ramped, because it cannot sample other, probably lower lying minima, from which its 
current state is separated by large (free energy) barriers. The complexity of the free energy 
landscape is usually greatly enhanced by the presence of disorder. It is well known [1-5], 
that disorder can lead to diverging barriers to relaxation and consequent nonequilibrium 
behavior and glassiness. 

(a) Avalanches: In some systems, collective behavior in the form of avalanches is found 
when the system is pushed by the driving field into a region of descending slope in the free 
energy surface. In experiments avalanches are often associated with crackling noises as in 
acoustic emission and Barkhausen noise [6-9]. There are other nonequilibrium systems where 
no such collective behavior is seen. Bending a copper bar for example causes a sluggish, 
creeping response due to the entanglement of dislocation lines. In contrast, wood snaps and 
crackles under stress due to "avalanches" of fiber breakings [10]. 

Although avalanches are collective events of processes happening on microscales, in many 
systems they can become monstrously large so that we - in spite of being large, slow creatures 



3 



- can actually perceive them directly without technical devices. This reminds one of the 
behavior observed near continuous phase transitions, where critical fluctuations do attain 
human length and time scales if a tunable parameter is close enough to its critical value. 
Correspondingly one might expect to find universal features when the sizes and times of the 
avalanches get large compared to microscopic scales. 

In fact, interesting questions concerning the distribution of avalanche sizes arise. Many 
experiments show power law distributions over several decades. For example, experi- 
ments measuring Barkhausen-pulses in an amorphous alloy, iron and alumel revealed several 
decades of power law scaling for the distribution of pulse areas, pulse durations and pulse- 
energies [11]. Similarly, Field, Witt, and Nori recorded superconductor vortex avalanches 
in Nb i7 %Ti 53 % in the Bean-state as the system was driven to the threshold of instability by 
the slow ramping of the external magnetic field. The avalanche sizes ranged from 50 to 10 7 
vortices. The corresponding distribution of avalanche sizes revealed about three decades of 
power law scaling. Numerous other systems show similar power law scaling behavior [11-15]. 

Why should there be avalanches of many sizes? Power laws suggest a scaling relationship 
between different length scales with universal exponents. 

There has been much recent progress studying avalanches near (continuous) depinning 
transitions. In these systems a single, preexisting interface or "rubber sheet" is pushed 
through a disordered medium by an external driving force. When the randomness is in 
some sense weak, the interface distorts elastically without breaking over a wide range of 
length scales [16]. The distortions occur in the form of avalanches on increasing size as 
the external driving force is raised to a critical threshold field at which the interface starts 
to slide ("depins"). Examples studied include charge density waves [17-22], weakly pinned 
Abrikosov flux lattices [23], single vortex-lines [24-26], preferentially wetting fluids invading 
porous media [27,28], a single advancing domain wall in weakly disordered magnets [28-30] 
and fluids advancing across dirty surfaces [31]. Their analytical description turns out to 
be rather involved, demanding functional renormalization groups (see appendix G). If the 
disorder is in some sense strong, the elastic interface can tear. It then responds much more 
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inhomogeneously and like a plastic or fluid, to the external driving force. This is the case for 
strongly pinned vortex lines in the mixed state of superconducting films [32,25], the invasion 
of nonwetting fluids into porous media [33,28] and nonlinear fluid flow across dirty surfaces 
[34] , and others. 

Many hysteretic systems exhibit a wide distribution of avalanche sizes without an under- 
lying depinning transition. They usually have many interacting advancing interfaces and in 
some cases also new interfaces created spontaneously in the bulk. We propose that the large 
range of observed avalanche sizes in these systems might be a manifestation of a nearby crit- 
ical point with both disorder and external magnetic field as tunable parameters. In the class 
of models, which we have been studying near the critical point, we find not only universal 
scaling behavior in the avalanche size distribution, but also in the shape of the associated 
hysteresis loops. 

(b) Hysteresis (response lags the force): Hysteresis is often observed at first order 
transitions, when the system goes out of equilibrium. It probably has been best studied in 
magnetic systems. The origin of magnetic hysteresis lies in the spontaneous magnetization 
of the microscopic Weiss-domains [35]. In the demagnetized state, the Weiss domains are 
irregularly oriented in different directions. The orientation of the magnetization of each 
domain is a function of the energy of the magnetic field and the elastic energy in the crystal 
(magnetostriction) and it is chosen such that the free energy is at a (local) minimum. Since 
there are many minima available, the specific choice depends on the history of the system. 
Thermal vibrations of the lattice are usually not sufficient to rotate the magnetization of 
a Weiss-domain into another preferred direction, since in most cases these directions are 
separated from each other by large energy barriers. As a weak magnetic field is applied, first 
those Weiss-domains which are most closely aligned with the field will grow at the expense of 
the others. At higher fields, up to saturation, entire domains will rotate in the direction of the 
field. The ultimate cause for hysteresis are the irreversible domain wall motion and domain 
rotation, which happen suddenly, as "avalanches" without further field increase, when the 
corresponding threshold fields are exceeded. The resulting jumps in the magnetization are 
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called Barkhausen jumps. Under a magnifying glass the magnetization curve looks like 
a staircase: the slope of the flat parts is due to the reversible part of the susceptibility, 
the step height is given by the irreversible avalanche- like changes of the magnetization. 
If the magnetization curve has the shape of a rectangle, the change of the magnetization 
happens in enormous, system sweeping avalanches, the so called large Barkhausen effect. 
Experimentally the Barkhausen jumps can be observed by magnetic induction or through 
the associated acoustic emission. 

Analogous effects are found in ferroelectric materials, where avalanches of flipping fer- 
roelectric domains can be observed in response to a changing external electric field [36,37]. 
Hysteresis curves with step-like noise are also found in elastic transformations, for exam- 
ple in athermal shape-memory alloys ramping temperature or stress. The noise is due to 
avalanches of regions transforming from martensite to austenite or vice versa [13]. Similar 
behavior has been observed for vortices moving in avalanches in type II superconductors as 
the external magnetic field is increased [12], for liquid helium leaving Nuclepore in avalanches 
as the chemical potential is reduced [15], and for some earthquake models [38-40,14]. In 
section II we discuss recent experiments performed in some of these systems. 

We have modeled the long wavelength, low frequency behavior of these systems using 
the nonequilibrium zero temperature random field Ising model (RFIM). Some of our results 
have been published previously [41,42]. In contrast to some other hysteresis models, like 
the Preisach model [43] and the Stoner-Wohlfarth model [7], where interactions between 
the individual hysteretic units (grains) are not included and collective behavior is not an 
issue, in the RFIM the intergrain coupling is the essential ingredient and cause for hysteresis 
and avalanche effects. Tuning the amount of disorder in the system we find a second order 
critical point with an associated diverging length scale, measuring the spatial extent of the 
avalanches of spin flips. 

(c) Avalanches in the RFIM: A power law distribution with avalanches of all sizes is 
seen only at the critical value of the disorder. However, our numerical simulations indicate 
that the critical region is remarkably large: almost three decades of power law scaling in 
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the avalanche size distribution remain when measured 40% away from the critical point. At 
2% away, we extrapolate seven decades of scaling. One reason for this large critical range is 
trivial: avalanche sizes are expressed in terms of volumes rather than lengths, so one decade 
of length scales translates to at least three decades of size (or more if the avalanches are not 
compact, i.e. if the Hausdorff dimension is less than three). Some experiments that revealed 
three decades of power law scaling have been interpreted as being spontaneously self-similar 
("self-organized critical") [12, 11, 44]. 1 Our model suggests that many of the samples might 
just have disorders within 40% of the critical value. Tuning the amount of disorder in these 
systems might reveal a plain old critical point rather than self organized criticality. 

(d) Hysteresis in the RFIM: At the critical disorder we also find a transition in the 
shape of the associated hysteresis loops: Systems with low disorder relative to the coupling 
strength, have rectangle-shaped hysteresis loops and a big (Barkhausen) discontinuity, while 
systems with large disorder relative to the coupling show smooth hysteresis loops without 
macroscopic jumps. At the critical disorder R c separating these two regimes, the size of the 
jump seen in the low disorder hysteresis loops shrinks to a point at a critical magnetic field 
H C (R C ), where the magnetization curve M(H) has infinite slope. The power law with which 
it approaches this point is universal. 

(e) Results: We have extracted the universal exponents near this transition point from 
a history dependent renormalization group (RG) description for the nonequilibrium zero 



1 The name "self-organized critical" is in fact an oxymoron: "Critical" means that you have to 
be just at the right place, "self-organized" indicates that the system does not need to be tuned 
anywhere special. While it is true that many "self organized critical" systems are regular critical 
points where the boundary conditions stabilize the system at the transition, the same is true of ice 
in a glass of water. The tools for studying the critical point are the RG methods in development 
over many years. We discuss here our doubts that the experiments are self-organized even in this 
limited sense. 
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temperature random field Ising model. The calculation turns out to be much simpler than 
for related depinning transitions [19,20,29-31,26] (see also appendix G). Above 6 dimensions 
the exponents are described by mean field theory. We expand the critical exponents around 
mean-field theory in 6 — e dimensions and discover a mapping to the perturbation expansion 
for the critical exponents in the pure equilibrium Ising model in two lower dimensions. The 
mapping does not, however, apply to the exponents governing the avalanche size distribution, 
which to our knowledge, have not yet been calculated directly in the depinning transitions. 
The simplicity of the RG calculation allowed us to develop a new method to calculate these 
avalanche exponents directly in the e-expansion, involving replicas of the system in a very 
physical way. We have used it to calculate the avalanche exponents to first order in e. We 
report under separate cover [45] extensive numerical simulations used to extract exponents in 
3, 4, and 5 dimensions (see section X). We find good agreement between the two approaches. 

This paper is organized as follows: In section II we discuss several experiments in mag- 
netic systems, shape memory alloys, porous media and superconductors that have close 
connections to the model studied here. The model is introduced in section III and a sum- 
mary of our results is given in section IV. In section V we pause for a moment and reflect 
upon our real motives. In section VI the RG description is set up using the Martin-Siggia- 
Rose formalism, and a description of the perturbative expansion and the results for the 
exponents to 0(e) is given in section VII. section VIII contains a discussion of the mapping 
of the expansion to the expansion for the thermal RFIM. We extract corrections to 0(e 5 ) 
for most of the exponents and show a comparison between the Borel resummation of the 
e-expansion and numerical results. In section IX a new method to calculate avalanche ex- 
ponents directly in an e-expansion is described and performed to 0(e). Finally, in section X 
we compare the results to our numerical simulation [45]. 

Some of the details of the mean-field calculation are given in appendix A. The expected 
tilting of the scaling axes in finite dimensions is discussed in appendix B. Details on the 
implementation of the history in the RG calculation are given in appendix C. Appendix D 
contains a description of the Borel summation of the results for rj and \ jv to 0(e 5 ) (which 
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is relevant also for the pure Ising model). The behavior near the infinite avalanche line in 
systems with less than critical randomness is discussed in appendix E. Appendix F renders 
details on the calculation of the avalanche exponents by the use of replicas. Related problems 
are finally discussed in appendix G. 

II. EXPERIMENTS 

In this section we will discuss several experiments that reveal scaling behavior which 
might be related to the critical point studied in this paper. The critical exponents found 
in real experiments do not necessarily have to be the same as in our model, since long- 
range interactions, different conserved quantities and other changes are likely to alter the 
universality class in some cases. We do propose however that the qualitative features, in 
particular the existence of an underlying plain old critical point with disorder and driving 
field as tunable parameters are likely to be the same as in our model (see also appendix G). 
A more detailled discussion of these and other Barkhausen experiments [91-96,99,97], and 
related experiments in non- magnetic avalanching systems (in shape memory alloys [13,46], 
superconductors [12,49], liquid helium in Nuclepore [15], and others), and a quantitative 
comparison with our theory will be given in a forthcoming publication [98,128,126]. 

A. Magnetic hysteresis loops for different annealing temperatures 

A beautiful, qualitative illustration of the crossover from smooth hysteresis loops at large 
disorder to hysteresis loops with macroscopic jumps at low disorder is shown in figure 1. 
The hysteresis loops were measured by Berger [47,48] for a 60 nm thick Gd film which had 
been grown onto a tungsten single crystal with a (110) surface orientation. The substrate 
as well as the film are highly purified (contaminants are less than 1/20 of a monolayer). Gd 
films prepared in this way exhibit a hep structure with the (0001) direction perpendicular to 
the surface. The substrate temperature during deposition was T = 350K, which results in 
smooth films with large atomically flat terraces, but also produces films with locally varying 
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strain and therefore with locally varying anisotropy. Subsequent annealing at higher tem- 
peratures improves the crystallographic order, which is accompanied by a strain relaxation. 
Thus, by varying the annealing temperature the authors are able to change the variation 
of the anisotropy defect density, which is somewhat analogous to the disorder parameter 
R in our model. Higher annealing temperatures correspond to lower values of R. If there 
is a second order critical point of the kind described in the introduction underlying the 
crossover from hysteresis loops with a jump to smooth hysteresis loops, it should be possible 
to extract a scaling form for the magnetization curves similar to the one given in eq. (9). 
(The annealing temperature minus some critical value would play the role of the tunable 
reduced disorder parameter r.) Under appropriate rescaling of the axes near the critical 
point the magnetization curves should all collapse one onto another. The necessary amount 
of rescaling as a function of distance from the critical point determines the (presumably 
universal) exponents [45]. Further measurements near the crossover to extract potential 
scaling behavior are currently being performed by Berger [47]. 

B. Barkhausen noise for different annealing temperatures 

Scaling behavior has also been recorded in the Barkhausen pulse duration and pulse 
area distribution in a related experiment [6]. The pulse area gives the total change in the 
magnetization due to the corresponding Barkhausen pulse. It is analogous to the avalanche 
size given by the number of spins participating in an avalanche in our model. It was found 
that the distribution of pulse areas integrated over the hysteresis loop of an 81% Ni-Fe 
wire (50 cm long, 1mm diameter) was well described by a power law up to a certain cutoff 
size (see figure 2). The cutoff appeared to be smaller at higher annealing temperatures. 
It would be interesting to see whether the cutoff takes a system-size dependent maximum 
value at a critical annealing temperature T c ann and decreases again at higher and lower 
annealing temperatures. This would be expected if varying the annealing temperature would 
correspond to tuning the system through a critical region with a diverging length scale at 
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T c ann . Near T c ann the Barkhausen pulse area distributions should then be described by a 
scaling form that would allow a scaling collapse of all distributions onto one single curve for 
appropriate stretching of the axes. Again, potentially universal critical exponents could be 
extracted from such a collapse. They would be predicted by our avalanche critical exponents 
if our model is in the same universality class. 

C. Barkhausen pulse size distributions at fixed disorder 

There are other experiments which revealed power law decays for Barkhausen pulse size 
distributions in various samples, as we have mentioned earlier [11,44]. To our knowledge the 
amount of disorder (or another parameter) was not varied in these experiments. The power 
law scaling over several decades found in these systems has in some cases been interpreted 
as a manifestation of self-organized criticality [11]. According to our simulations however, 
three decades of scaling occur when the disorder is as far as 40% away from the critical 
value. Tuning the amount of disorder in the system (for example by annealing the sample 
as in the previous two experiments, or by introducing random strain fields) might lead to a 
larger (or smaller) cutoff in the power law pulse size distribution. It seems rather plausible 
that the observed scaling behavior would be due to a plain old critical point rather than 
self-organized criticality. 

D. Remarks 

The first experiment that showed the crossover in the shape of the magnetic hysteresis 
loops was performed in an effectively two dimensional system, while the experiments on 
Barkhausen noise used effectively three dimensional systems. Interestingly, two might be 
the lower critical dimension of the transition which we are studying in this paper. 2 These 



2 In one dimension there will still be a crossover from hysteresis loops with a macroscopic jump 
to smooth hysteresis loops for a bounded distribution of random fields [51], however the potential 
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conjectures are currently being tested with our numerical simulations [45]. 

Real experiments may involve long-range fields which may in principle alter the uni- 
versality class. Different kinds of disorder, such as correlated disorder rather than point 
disorder and random anisotropies or random bonds rather than random fields, may also be 
present. Furthermore, the symmetries can be changed if there are more than two available 
("spin") states at each site in the lattice. In appendix G we discuss related models, and 
which of these changes are expected to change the universality class relative to our model. 

What is the moral of this story? Power laws with cutoffs are anything but sufficient 
evidence for self-organized criticality. In the systems discussed here, they are more likely 
due to a nearby plain old critical point with disorder as a tunable parameter and a large 
critical region, than spontaneous self-organization towards a self-similar state. For related 
future experiments and analysis one would recommend the search for tunable parameters 
other than system size that allow to change the cutoff in the power law distributions. For the 
analysis of our simulation results, scaling collapses and other techniques from equilibrium 
critical phenomena turned out to be very useful for extracting critical exponents [45]. 

III. THE MODEL 

As we have explained in the introduction, the goal is to describe the long-wavelength 
behavior of hysteretic systems with noise due to microscopic avalanches triggered by the 
external driving field. We will focus in particular on the scaling regime, where collective 
behavior is observed on many length scales. Conventional hysteresis models like the Preisach 
model [43], which do not take into account interaction between the smallest hysteretic units 
(grains), would not be suitable for this purpose. The Preisach model could only be used 
to fit a certain measured distribution of avalanche sizes — the power law scaling would be 



scaling behavior found near the transition will not be universal, but rather depend on the exact 
shape of the tails of the distribution of random fields. 
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the input determining free parameters of the model rather than the output with universal 
predictive power. 

The key ingredient is interaction. 3 As is well known from equilibrium phenomena, behav- 
ior on long length scales can often be well described by simple microscopic models that only 
need describe a few basic properties correctly, such as symmetries, interaction range and 
effective dimensions. This notion has been successfully applied in particular to equilibrium 
magnetic systems: the scaling behavior found in some pure anisotropic ferromagnets near 
the Curie-temperature is mimicked reliably by the regular Ising model [53,54]. At each site 
% in a simple cubic lattice there is a variable Sj, in this context called a spin, which can take 
two different values, Sj = +1 or s$ = — 1 [55]. (This corresponds to a real magnet where a 
crystal anisotropy prefers the magnetic moments (spins) to point along a certain easy axis.) 
Each spin interacts with its nearest neighbors on the lattice through an exchange interac- 
tion, Jij = J/z, which favors parallel alignment, z is the coordination number of the lattice 
and J is a positive constant. (For the behavior on long length scales the exact range of the 
microscopic interaction is irrelevant, so long as it is finite.) One can write the Hamiltonian 
as 

H = -^JijSiSj - H^Si, (1) 

ij i 

where it is understood that the sum runs over nearest neighbor pairs of spins on sites i and 
j. H is a homogeneous external magnetic field. In two and higher dimensions this model 



3 For systems exhibiting return-point-memory (also called "subloop-closure" or "wiping-out" prop- 
erty) [41,52,43] there is in fact a well established method to verify whether interactions play an 
important role, and whether the Preisach model is applicable at all. It involves testing for subloop- 
congruency. For further details we refer the reader to the literature [52,43,15]. This test has 
been conducted for the experiment on liquid He in Nuclepore [15] revealing the importance of 
interactions between the pores, and the failure of the conventional Preisach model to describe the 
hysteresis curve. 
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exhibits an equilibrium ferromagnetic state at temperatures T < T c , where T c is the Curie 
temperature. Figure 3 shows the corresponding equilibrium magnetization curve at zero 
temperature: All spins are pointing up at positive external magnetic fields, and all spins are 
pointing down at negative external magnetic fields. At if = the curve is discontinuous. 

What would the magnetization curve look like for the same model, but far from equi- 
librium, as is the case for most real magnets? The answer is shown in figure 4. We have 
imposed a certain local dynamics, assuming that each spin s« will flip only when the total 
effective field at its site, given by 

hf^-^JijSj-H, (2) 

3 

changes sign. We find that the resulting magnetization curve becomes history dependent. 
The system will typically be in some metastable state rather than the ground state. The 
upper branch of the hysteresis curve in figure 4 corresponds to the case where we have 
monotonically and adiabatically lowered the external magnetic field, starting from if = +oo, 
where all spins were pointing up. At the coercive field H l c = —2dJij = — J all spins flip in 
a single system spanning event or "avalanche". Similarly, for increasing external magnetic 
field, they all flip at if" = 2dJ„ = +J. It becomes clear that the underlying cause for 
hysteresis in this model is the interaction between the spins. 

So far, however, an essential feature of real materials is missing: there is no account for 
dirt. Usually there will be inhomogeneities and disorder in the form of defects, grain bound- 
aries, impurities, leading to random crystal anisotropies, and varying interaction strengths 
in the system. Consequently not all spins will flip at the same value of the external magnetic 
field. Instead, they will flip in avalanches of various sizes that can be broken up or stopped 
by strongly "pinned" spins or clusters of previously flipped spins. 

If the disorder in the system is small, the picture will not deviate dramatically from the 
pure case. One would expect only a few small precursors to the macroscopic avalanche of 
figure 4. If however, the disorder is large compared to the coupling strength in the system, 
one might expect no system sweeping avalanche at all, but only small clusters of spins 
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flipping over a broad range of the external magnetic field. 

A simple way to implement a certain kind of uncorrelated, quenched disorder is by 
introducing uncorrelated random fields into the model. (Other kinds of disorder are discussed 
in appendix G.) The energy function is replaced by 

H = ~J2 J ij S i S j - J2( H + fi) S i- ( 3 ) 

ij i 

The local dynamics remains unchanged, except for a modification in the expression for the 
total effective field at site i, which now also has to take into account the random field ff 



hf = -E J a 8 i- H -fi- ( 4 ) 

3 



We assume a Gaussian distribution p(fi) of standard deviation R for the fields which is 
centered at fi — 0: 

As we will show, the critical exponents do not depend on the exact shape of the distribution 
of random fields. To pick a Gaussian is a standard choice, which (due to the central limit 
theorem [56]) is also more likely to be found in some real experiments than, for example, 
rectangular distributions [57]. 

In magnets the random fields might model frozen-in magnetic clusters with net magnetic 
moments that remain fixed even if the surrounding spins change their orientation. In con- 
trast to random anisotropies they break time reversal invariance by coupling to the order 
parameter (rather than its square). In shape memory alloys, ramping temperature, the ran- 
dom fields can be thought of as concentration fluctuations that prefer martensite over the 
austenite phase [58]. In the martensitic phase, ramping stress, they model strain fields that 
prefer one martensitic variant over another [58]. 

In appendix G we discuss related systems with different kinds of disorder, symmetries, 
interactions and dynamics, and the possible effects of such changes on the associated long- 
wavelength behavior and critical properties at transitions analogous to the one studied here. 
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IV. RESULTS 



It is relatively easy to solve the nonequilibrium model of eq. (3) in the mean-field ap- 
proximation where every spin interacts equally strongly with every other spin in the system. 
The coupling is of size = J/N, where N is the total number of spins 4 (i.e. all spins act 
as nearest neighbors). The Hamiltonian then takes the form 

H = -^2(JM + H + f i )s i , (6) 

i 

Just as in the Curie- Weiss mean field theory for the Ising model, the interaction of a spin 
with its neighbors is replaced by its interaction with the magnetization of the system. 

It turns out that the mean-field theory already reflects most of the essential features of 
the long-length scale behavior of the system in finite dimensions: Sweeping the external field 
through zero, the model exhibits hysteresis. As disorder is added, one finds a continuous 
transition where the jump in the magnetization (corresponding to an infinite avalanche) 
decreases to zero. At this transition power law distributions of noise (avalanches) and 
universal behavior are observed. 

As we will show later in an RG description of the model, the critical exponents describing 
the scaling behavior near the critical point are correctly given by mean-field theory for 
systems in 6 and higher spatial dimensions. The RG allows us to calculate their values 
in (6 — e) dimensions in a power series expansion in e > around their mean-field values 
at e = 0. In the following we briefly present the results from mean-field theory, from the 
e-expansion and from numerical simulations in 3 dimensions. More details will be given in 
later sections. 



4 We warn the reader that in all plots of numerical simulation results in finite dimensions J will 
denote the strength of the nearest neighbor coupling Jjj in the simulated crystal, while in the 
analytic calculation and in mean-field theory it denotes J2j Jij-> i- e - the two definitions differ by 
the coordination number of the lattice. 
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A. Results on the magnetization curve 



Figure 5 shows the hysteresis curve in mean field theory at various values of the disorder 
R < R c = yj (2/n)J, R = R c , and R > R c . For R < R c , where the coupling is important 
relative to the amount of disorder in the system, the hysteresis curve displays a jump due 
to an infinite avalanche of spin flips, which spans the system. Close to R c the size of the 
jump scales as AM ~ r 13 , with r = (R c — R)/R, and (3 = 1/2 in mean field theory. Using a 
mapping to the pure Ising model, we find in 6 — e dimensions [59] 

/? = 1/2 — e/6 + 0.00617685e 2 - 0.035198e 3 + 0.0795387e 4 

-0.246111e 5 + O(e 6 ). (7) 

At R — R c the magnetization curve scales as M — M(H C (R C )) ~ h^^ 5 \ where h = 
H—H C (R C ) and H C (R C ) is the (nonuniversal) magnetic field value at which the magnetization 
curve has infinite slope. In this mean field theory H C (R C ) = 0, and M(H C (R C )) = 0, and 
(35 — 3/2. In 6 — e dimensions [59] 

(35 = 3/2 + 0.0833454e 2 - 0.0841566e 3 + 0.223194e 4 - 0.69259e 5 + 0(e 6 ) . (8) 

Numerical simulations in 3 dimensions yield (3 = 0.036 ± 0.036 and (35 = 1.81 ± 0.36 [45]. 

For R > R c the disorder can be considered more important than the coupling. Conse- 
quently there are no system spanning avalanches (for infinite system size) and the magneti- 
zation curve is smooth. 

Note that the hard spin mean- field theory does not show any hysteresis for R > R c . 
This is only an artifact of its particularly simple structure and not a universal feature. For 
example, the analogous soft spin model, which is introduced for the RG description in section 
VI, has the same exponents in mean-field theory, but shows hysteresis at all disorders R, 
even for R > R c as seen in figure 6. 

Close to R c and H C (R C ) the magnetization curve is described by a scaling form: 

M - M(H C (R C )) = m(r, h) ~ r f3 M±(h/r l3S ) (9) 
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where A4± is a universal scaling function (± refers to the sign of r). It is computed in 
mean-field theory in appendix A. Corrections 5 to the mean-field equation of state in 6 — e 
dimensions are calculated to 0(e) in section VII B 4. Results to 0(e 2 ) are quoted in section 
VIII. 

B. Results on the mean-field phase diagram 

Figure 7 shows the phase diagram for the lower branch of the hysteresis curve as obtained 
from the simple hard spin mean-field theory, defined through eq. (6). The bold line with 
the critical endpoint (R c , H C (R C )) indicates the function H^(R) for the onset of the infinite 
avalanche for the history of an increasing external magnetic field. The dashed line describes 
H l c (R) for a decreasing external magnetic field. The three dotted vertical lines marked (a), 
(b), and (c) describe the paths in parameter space which lead to the corresponding hysteresis 
loops shown in figure 5. Figure 8 shows the corresponding phase diagram for the soft-spin 
model. As before, the three dotted vertical lines marked (a), (b) and (c) indicate the paths 
through parameter space associated with the three hysteresis loops shown in figure 6. Note 
that in the soft-spin model as well as in simulations in finite dimensions the two infinite 
avalanche lines H^{R) and H l c (R) do not touch at R c — this is another way of saying that 
in these cases there is hysteresis at R = R c , and, because of continuity also at R > R c . 

C. Results on scaling near the onset H C (R) of the infinite avalanche line (R < R c ) 

The mean-field magnetization curve scales near the onset of the infinite avalanche as 



5 The scaling form in finite dimensions may depend not on r and h, but on rotated variables r' 
and h', which are linear combinations of r and h. This applies to all scaling relations derived in 
the infinite range model. The amount by which the scaling axes r' = and h' = are turned 
relative to r = and h = is a nonuniversal quantity and has no effect on the critical exponents 
(see appendix B). 
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(M - M C {H C {R))) ~ (if - H C (R))< 



(10) 



with ( — 1/2. (In the following H C (R) always stands for H™(R) for the history of an 
increasing external magnetic field, and for H l c (R) for a decreasing external magnetic field.) 

Curiously we do not observe this scaling behavior in numerical simulations with short 
range interactions in 2, 3, 4 and 5 dimensions. Indeed, the RG description suggests that 
the onset of the infinite avalanche would be an abrupt ("first order" type) transition for all 
dimensions d < 8 (see appendix E), and a continuous transition for d > 8. We [45] have 
performed initial numerical simulations in 7 and 9 dimensions for system sizes 7 7 and 5 9 
at less than critical disorders. The simulation results do in fact seem to confirm the RG 
prediction [45]. In the following we will mostly focus on the critical endpoint at (R c , H C (R C )), 
where the mean-field scaling behavior is expected to persist in finite dimensions with slightly 
changed critical exponents. 



Figure 9 shows the magnetization curves from simulations of two 3 dimensional systems 
with only 5 3 spins. The curves are not smooth. They display steps of various sizes. Each step 
in the magnetization curve corresponds to an avalanche of spin flips during which the external 
magnetic field is kept constant. Figure 10 shows histograms D(S,r) of all avalanche sizes 
S observed in mean-field systems at various disorders r when sweeping through the entire 
hysteresis loop. For small r the distribution roughly follows a power law D(S, r) ~ g-ir+eps) 
up to a certain cutoff size S max ~ |r| _1 / CT which scales to infinity as r is taken to zero. 

In appendix A we derive a scaling form for the avalanche size distribution for systems 
near the critical point: Let D(S,r, h) denote the probability to find an avalanche of size S 
in a system with disorder r at magnetic field h upon an infinitesimal increase of the external 
magnetic field. For large S one finds 



D. Results on avalanches 



D(S,r, h) ~ l/S T V±(Sr 



(11) 
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The scaling form for D(S,r) of the histograms in figure 10 is obtained by integrating 
D(S,r, h) over the external magnetic field. 

In mean field theory we find a = 1/2 and t = 3/2. In 6 — e dimensions we obtain from 
the RG calculation a — 1/2 — e/12 + 0(e 2 ) and r = 3/2 + 0(e 2 ). Numerical simulations in 
3 dimensions [45] render a = 0.238 ± 0.017 and r = 1.60 ± 0.08 [45]. 

E. Results on correlations near (R C ,H C (R C )) 

With the mean field approximation we have lost all information about length scales in the 
system. The RG description, which involves a coarse graining transformation to longer and 
longer length scales, provides a natural means to extract scaling forms for various correlation 
functions and the correlation length of the system. 

1. Avalanche correlations 

The avalanche correlation function G(x,r, h) measures the probability for the configura- 
tion of random fields in the system to be such that a flipping spin will trigger another at 
relative distance x through an avalanche of spin flips. Close to the critical point and for 
large x the function G(x, r, h) scales as 

G(x,r,h) ~ l/x d - 2 ^g ± {x/i{r,h)) , (12) 

where 77 is called "anomalous dimension" and Q± is a universal scaling function. The correla- 
tion length £(r, h) is the important (macroscopic) length scale of the system. At the critical 
point, where it diverges, the correlation function G(x, 0, 0) decays algebraically — there will 
be avalanches on all length scales. Close to the critical point the correlation length scales as 

ar,h)~ r -»y ± (h/rP s ), (13) 

where y± is the corresponding scaling function. From the e-expansion one obtains [59] 

\/v = 2 - e/3 - 0.1173e 2 + 0.1245e 3 - 0.307e 4 + 0.951e 5 + 0(e 6 ) , (14) 
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and 

r] = 0.0185185e 2 + 0.01869e 3 - 0.00832876e 4 + 0.02566e 5 + 0(e 6 ) . (15) 
The numerical values in 3 dimensions are \ jv — 0.704 ± 0.085 and rj = 0.79 ± 0.29. 

2. Spin-spin ("cluster") correlations 

There is another correlation function which measures correlations in the fluctuations of 
the spin orientation at different sites. It is related to the probability that two spins and 
Sj at two different sites i and j, that are distanced by x, have the same value [53]. It is 
defined as 

C(x,r,h) = (( Si - {si) f )(sj - (sj)f)) f , (16) 

where () / indicates the average over the random fields. From the RG description we find 
that for large x it has the scaling form 

C(x,r,h) ~x-( d - 4+ ^C±(x/ar,h)), (17) 

where £(r, h) scales as given in eq. (13) and C± is a universal scaling function. At the critical 
point C(x, 0, 0) decays algebraically — there will be clusters of equally oriented spins on all 
length scales. The e-expansion renders [59] 

fj = 0.0185185e 2 + 0.01869e 3 - 0.00832876e 4 + 0.02566e 5 + 0(e 6 ) . (18) 

which is in fact the same perturbation expansion as for rj to all orders in e. (The two 
exponents do not have to be equal beyond perturbation theory, see also [126,129]). 

F. Results on avalanche durations 

Avalanches take a certain amount of time to spread, because the spins are flipping 
sequentially. The further the avalanche spreads, the longer it takes till its completion. The 
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RG treatment suggests that there is a scaling relation between the duration T of an avalanche 
and its linear extent / 



T{1) ~ r 



(19) 



with z = 2 + 2r] to 0(e 3 ) [60], i.e. 



z = 2 + 0.037037e 2 + 0.03738e 3 + 0(e 4 ) . 



(20) 



Our numerical result in 3 dimensions is z — 1.7 ± 0.3. 6 



The fractal dimension for the biggest avalanches S, 



max 



IS 



d fractal = , 



(21) 



so that the time for the biggest finite avalanches scales as T(S max ) ~ S\ 



<auz 
max ' 



G. Results on the area of the hysteresis loop 



In some analogy to the free energy density in equilibrium systems, one can extract the 
scaling of the area of the hysteresis loop for this system near the critical endpoint. (This is 
the energy dissipated in the loop per unit volume.) From the fact that the singular part of 
the magnetization curve scales as m(h,r) ~ r f3 A4±(h/r 135 ) (see eq. (9)) we conjecture that 
the singular part of the area would scale as A sing ~ / m(h, r)dh ~ r 2 ~ a with 2 — a = (3 + (58. 
(The scaling form for the total area A tot will also have an analytical piece: A tot = c + Cir™ + 



6 While we expect the 6 — e results for the static exponents /3, S, u, rj, fj, r, a, etc. to agree 
with our hard-spin simulation results close to 6 dimensions, this is not necessarily so for the 
dynamical exponent z. There are precedences for the dynamics being sensitive to the exact shape 
of the potential, sometimes only in mean-field theory [17,19,20] (charge density waves in a smooth- 
potential versus a linear cusp potential), and sometimes even in the e-expansion [20,127] (charge 
density waves in a sawtooth potential). 
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c 2 r( n+1 ) + ... + A sing ; for any data analysis one needs to keep all terms with n < 2 — a.) In 
mean field theory a = 0. Numerical and analytical results can be derived from the results 
for j3 and (38 quoted earlier. 

H. Results on the number of system-spanning avalanches at the critical disorder 

R = R c 

In percolation in any dimensions less than 6, there is at most one infinite cluster present 
at any value of the concentration parameter p, in particular also at its critical value p c 
[61]. In contrast, in our system at the critical point R = R c the number iVoo of "infinite 
avalanches" found during one sweep through the hysteresis loop, diverges with system size 
as iV^ ~ L e in all dimensions d > 2. 

The e-expansion for our system yields 

9v = 1/2 - e/6 + 0(e 2 ) (22) 

Numerical simulations [45] show clearly that 9 > in 4 and 5 dimensions. In three di- 
mensions one finds Qv = 0.021 ± 0.021 and 9 = 0.015 ± 0.015. (For more details on 9 see 
[126,129].) 

I. List of exponent relations 

In the following sections we list various exponent relations, for which we give detailled 
arguments in references [126,129]. 

1. Exponent equalities 
The exponents introduced above are related by the following exponent equalities: 

P-P5 = {t-2)/<j if r<2, (23) 
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(24) 



P = z(d-4 + fj), (25) 

and 

5= (d-2ri + fj)/(d-4 + fj). (26) 

(The latter three equations are not independent and are also valid in the equilibrium random 
field Ising model [62,125,101]). 7 

2. Incorrect exponent equalities 

a. Breakdown of hypers caling In our system there are two different violations of hyper- 
scaling. 

1. In references [126,129], we show that the connectivity hyperscaling relation 1/a — 
dv — f3 from percolation is violated in our system. There is a new exponent 9 defined by 
1/a = (d-9)v-f3 with 9v = 1/2 - e/6 + 0(e 2 ) and 9v = 0.021 ±0.021 in three dimensions 
[45]. 9 is related to the number of system spanning avalanches observed during a sweep 
through the hysteresis loop. 

2. As we will discuss in section VIII there is a mapping of the perturbation theory for 
our problem to that of the equilibrium random field Ising model to all orders in e. From that 
mapping we deduce the breakdown of an infamous ( "energy" )-hyperscaling relation, which 
has caused much controversy in the case of the equilibrium random-field Ising model [62] 

P + p6=(d-9)v, (27) 



7 Also, using these relations one finds that the inequality v j fib > 2/d (which applies in our model 
[129,126,104]) goes over into the Schwartz-Soffer inequality fj < 2n that has been derived for the 
corresponding equilibrium model [63]. 
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with a new exponent 9. In [126,129] we discuss the relation of the exponent 9 to the energy 
output of the avalanches. The e-expansion yields 9 = 2 to all orders in e. Non-perturbative 
corrections are expected to lead to deviations of 9 from 2 as the dimension is lowered. The 
same is true as is the case in the equilibrium RFIM [126,129]. The numerical result in three 
dimensions is 9 = 1.5 ± 0.5 [45]. (In the three-dimensional Ising model it is 9 eq = 1.5 ± 0.4 
[62,100].) 

b. Breakdown of perturbative exponent equalities There is another strictly perturbative 
exponent equality, which is also obtained from the perturbative mapping to the random-field 
Ising-model [126,129], 

77 = 77. (28) 
It, too, is expected to be violated by non-perturbative corrections below 6 dimensions. 

3. Exponent inequalities 

In references [126,129] we give arguments for the following two exponent-inequalities 8 : 

u/(35>2/d, (29) 

which is formally equivalent to the "Schwartz-Soffer" inequality, fj < 2i], first derived for the 
equilibrium random field Ising model [63], and 

v>2/d, (30) 

which is a weaker bound than eq. (29) so long as (55 > 1, as appears to be the case both 
theoretically and numerically at least for d > 3. 



8 From the normalization of the avalanche size distribution D(s,r,h) (see eq. (11)) follows that 
r > 1. 
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J. Results on the upper critical dimension of the critical endpoint (R C: H C (R C )) 



The consistency of the mean-field theory exponents for d > 6 can be shown by a Harris 
criterion type of argument [53], which also leads to eq. (29) [128,126]. Approaching the 
critical point along the r' — line, one finds a well defined transition point only if the 
fluctuations 5h' in the critical field H c due to fluctuations in the random fields are always 
small compared to to the distance hi from the critical point, i.e. Sh'/h' « 1 as h! — > 0. 
With 6ti ~ r d/2 and f ~ {r')- v f±{h'/{r'Y 5 ) ~ {h')- u/m at r' = 0, one obtains Sh'/h! ~ 
£-d/2^-f3S/u << 1; or ^5 > 2 /d. For v = 1/2 and (35 = 3/2 this is only fulfilled if d > 6, 
ie. d = 6 is the upper critical dimension. 

V. WHY AN e-EXPANSION? 

Perturbation theory has proven an invaluable tool for practical calculations in many 
branches of physics. The asymptotic expansion in powers of the electronic charge in elec- 
trodynamics is an example where perturbation theory not only gives extremely accurate 
results but also provides a qualitative insight in the underlying physical processes such as 
absorption and emission of photons. Critical phenomena in general are among the cases 
where perturbation theory cannot be applied (at least not below a certain "upper critical 
dimension" d c ). The fluctuations in the order parameter near the critical point become too 
large to extract useful information from a perturbation expansion in some physical coupling 
constant. 

The RG has been developed specifically for such cases. It is a means to extend the 
derivation of critical scaling forms from mean-field theory to finite dimensions, and to obtain 
information about the effective behavior of the system on long length scales. The Wilson- 
Fisher momentum shell renormalization group is an iterative coarse graining transformation. 
In each step the shortest wavelength degrees of freedom are integrated out, leading to an 
effective action for the modes on longer length scales. If a theory is renormalizable, the 
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coarse grained action will be in the same form as the original action, but with rescaled 
parameters. 

Under coarse graining many of these parameters will flow to a certain fixed point, which is 
independent of their original microscopic value. The properties of the system on long length 
scales only depend on the fixed point. This is the key to universality: different microscopic 
systems flow to the same fixed point and are therefore described by the same effective action 
on long length scales. They will consequently show the same critical behavior. 

Second order critical points are fixed points of the coarse graining transformation. They 
are characterized by a diverging correlation length and consequent scale invariance. There 
are fluctuations (such as clusters, avalanches, etc.) on all length scales. In the same sense 
power laws are scale free functions. 9 It is therefore not surprising that systems near their 
critical point are described by power laws. One of the triumphs of the RG is to show from 
first principles that near the critical point the interesting long-wavelength properties are 
given by homogeneous functions 10 with respect to a change of length scale in the system. 
This observation leads to Widom scaling forms for the various macroscopic quantities, some 
of which we already obtained from simple expansions in mean field theory. 

The RG thus provides a formal justification of the scaling ansatz used in the data anal- 
ysis and an explanation for the universality of the critical exponents. It is the ultimate 



9 When a function f(x) ~ is measured over three pairs of octaves, say over the intervals [1,4], 
[10,40], and [100,400], the ratio of the largest to the smallest value is always 4l a l, so the three 
graphs of f(x) can be superimposed by simple change of scale. In this sense, power laws are scale 
invariant [64]. 

10 A function f(xi, ■■■,x n ) is homogeneous of degree D in the variables x±,---,x n , if on multiplying 
each Xi by an arbitrary factor b the value of / is multiplied by b D , i.e. f(bx±, • • -, bx n ) = b D f(x±, • • 
■,x n ) [64]. Here we use the definition f(b Xl xi,- ■ -,b x "x n ) = b y f(x\,- ■ -,x n ) for Ai,- • -,A n real 
constants. 
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justification for the attempt to extract useful predictions about real complex materials from 
extremely simple caricatures of the microscopic physics. 

In particular it has been used to derive an expansion for the critical exponents around 
their mean-field values in powers of the dimensional parameter e = d c — d, where d is the 
dimension of the system. Note that the e-expansion is an asymptotic expansion in terms of 
a quite unphysical parameter. Nevertheless it has proven very successful for mathematical 
extrapolations. In this paper we shall apply its basic ideas to our problem and refer the 
reader for further details to excellent reviews in the existing literature [65,66,53,64,67-70]. 

The calculation turns out to be interesting in its own right. In contrast to RG treat- 
ments of equilibrium critical phenomena, a calculation for our hysteresis problem has to 
take into account the entire history of the system. It reveals formal similarities to related 
single interface depinning transitions [19,20,29-31,26]. Although our problem deals with the 
seemingly more complex case of many advancing interfaces or domain walls, the calculation 
turns out to be rather simple, much simpler in fact than in the single interface depinning 
problem. More details are given in appendix G. The techniques employed here are likely to 
be applicable to other nonequilibrium systems as well. 

VI. ANALYTICAL DESCRIPTION 

In equilibrium systems one defines a partition function as the sum over the thermal 
weights or probabilities of all possible states. From this partition function all ensemble 
averaged correlation functions can be obtained. It is also usually the quantity used to 
calculate the critical exponents in an RG treatment. Our system is at zero temperature 
and far from equilibrium. For a given configuration of random fields, the system will follow 
a deterministic path through the space of spin microstates as the external magnetic field 
is raised adiabatically. Systems with different configurations of random fields will follow 
different paths. If we assign a 5-function weight to the correct path for each configuration 
and then average over the distribution of random fields, we obtain a probability distribution 
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for the possible paths of the system, which is the analogue of the probability distributions 
for the possible states in equilibrium systems. The analogue of ensemble averaging for 
equilibrium systems, is random field averaging in our system. The sum over all possible paths 
weighted by their corresponding probability will play the role of a partition function for our 
nonequilibrium system when we set up a Wilson-Fisher momentum shell renormalization 
group transformation to calculate the critical exponents. 

How can one formally describe the path that a system takes for a specific configuration 
of random fields as the external magnetic field is increased from — oo to +00? A convenient 
way is to introduce a time t into the otherwise adiabatic problem via H(t) — H + Qt. 
H is the magnetic field at time t — 0. Q > is the sweeping rate for a monotonically 
increasing external magnetic field. The idea is to write down an equation of motion for 
each spin such that the resulting set of coupled differential equations has a unique solution 
which corresponds to the correct path the system takes for a given history. By taking fl to 
zero in the end one obtains the adiabatic or "static" limit, in which we are interested. For 
convenience we introduce soft spins that can take values ranging from —00 to +00. Later 
on this will allow us to replace traces over all possible spin configurations by path integrals 
over the range of definition of the spins. We assume that each spin is moving in a double 
well potential V(si) with minima at the "discrete" spin values Sj = ±1: 



k/2 (si + l) 2 for s < 
' 1 ' (31) 

k/2 ( Si - l) 2 for s > 



To guarantee that the system takes a finite magnetization at any magnetic field, one needs 
k > and k/J > 1. The Hamiltonian of the soft spin model is then given by: 

W = - E J H s i s i - E(M + Hs i ' V ( s i)) ■ ( 32 ) 

ij i 

All terms in the Hamiltonian are as before, except for the additional V(sj) term. A spin flip 
in this model corresponds to a spin moving from the lower to the upper potential well. 11 



We believe the calculation could just as well have been performed for discrete spins, maybe 
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We impose purely relaxational dynamics, given by 

(l/T )d tSl (t) = -8H/8 Si (t) . (33) 

r is a "friction constant" . 

This model shows qualitatively similar behavior to real magnets: As the external mag- 
netic field is ramped, we observe spin flips, which correspond to irreversible domain wall 
motions. The linear relaxation between the spin flips corresponds to the reversible domain 
wall motion, which we have described in the introduction. 

The soft spin mean-field theory, where every spin interacts equally with every other spin 
yields the same static critical exponents as we have obtained earlier for the hard spin model. 
We have also checked that replacing the linear cusp potential by the more common, smooth 
s 4 double well potential does not change the static mean field exponents. 12 

A. Formalism 

We use the formalism introduced by Martin Siggia and Rose [73] to treat dynamical 
critical phenomena, which is similar to the Bausch-Janssen- Wagner method [74]. One 
defines the generating functional Z for the dynamical problem as an integral over a product 
of 5-functions (one for each spin), each of which imposes the equation of motion at all times 
on its particular spin [19]: 



by using the ideas of Lubensky et al. [71]. At the time it seemed easier to use the established 
formalism for continuous spins. However, in mean-field theory the formulas look very similar to 
the hard-spin model. The new parameter k from the double well potential never seems to come 
into play for any of the universal properties of the system. 

12 The form of the potential does change the (dynamical) mean-field exponents in the CDW de- 
pinning transition at zero temperature [19]. In finite dimensions it changes some properties that 
are associated with the thermal rounding of the CDW transition [72] . 
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(34) 



[ds] symbolizes the path integral over all spins in the lattice at all times, and J[s] is the 
necessary Jacobian, which fixes the measure of the integrations over the such that the 
integral over each delta-function yields 1 [19]. One can show that J[s] merely cancels the 
equal time response functions [75, 19]. 13 

In order to write Z in an exponential form in analogy to the partition function in equi- 
librium problems, we express the (^-functions in their Fourier-representation, introducing an 
unphysical auxiliary field Sj(t): 



Absorbing any constants into J[s\, this yields for the (not yet random-field averaged) gen- 
erating functional (in continuous time): 



with the action 



13 To that end, one chooses the following regularization (when discretizing in time) : Let t = ne with 
e a small number taken to zero later when n is taken to infinity, such that their product remains 
fixed. Then dtSi(t) becomes (sj(rae) — Sj((n — l)e))/e and we integrate over dsj(ne). Because of 
the analyticity we are free to take the rest of the argument of the (5-function at the lower value 
Sj((n — l)e), or the average value (sj(ne) + Sj((n — l)e))/2, or the upper value Sj(rae). If we choose 
the first possibility, one finds that the Jacobian will be only a constant, since each argument only 
depends on Sj(ne), and is independent of any other Sj{ne) with i ^ j. This corresponds to allowing 
a force at time (n — l)e to have an effect only after some time e, i.e. equal times response functions 
are manifestly zero. 





(36) 



or 




(37) 
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W = iJ2 I dtSjit^dtSj^/Fo + 6H/6sj(t)) (38) 

3 

= i ]T [ dts^dts^/To - ]T Jjisi ~H-f j + SV/Ssj) . 
j J I 

We can express correlation and response functions of Sj(t) as path integrals in terms of 
W, because solely the unique deterministic path of the system for the given configuration 
of random fields makes a nonzero contribution to the path integral over [ds] in eq. (36). For 
example the value of spin Sj at time t' is given by 

Sj (t') = Z- 1 j j [ds'Wds'Wmt) exp(W) . (39) 

Similarly, correlation functions are given by 

Sj (t')s k (t") = Z' 1 J J [ds'\[ds'\J[s%{t')J k {t") exp(W0 . (40) 

To calculate the response of Sj at time t' to a perturbative field Jek(t', t") switched on at site 
k at time t", we add the perturbation to the magnetic field at site k, such that the action 
becomes 

W e = iJ2 I dts 3 (t){d tSj (t)/T Q - ]T J jlSl -H-fj + 5V/5 Sj ) (41) 
j+k J i 

+i [ dts k (t)(d t s k (t)/T -J2JkiSi- H - f k + SV/Ss k - Je k (t,t")). 
J i 

Taking the derivative with respect to Je k and the limit e k — > afterwards 14 one obtains 

5 Sj (t')/5e k (t") = {-%)Z' X J J [ds}[ds}J[ S } Sj (t')s k (t") exp(W) , (42) 

so s acts as a "response field". Henceforth we shall suppress J\s\, keeping in mind that its 
only effect is to cancel equal time response functions (see previous footnote). 

Since Z — 1 independent of the random fields, we could have left out the Z~ l factors 
in eqs. (39), (40), and (42). This greatly facilitates averaging over the random fields: The 



1 An exact definition of the functional derivative which is consistent with the history of a mono- 
tonically increasing magnetic field is given in appendix C. 
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average response and correlation functions are generated by averaging Z directly over the 
random fields. Unlike in equilibrium problems with quenched randomness it is not necessary 
to calculate the (more complicated) average of InZ. One obtains for the random field 
averaged correlation functions 

<*;(*>*(*")>/ = / / [ds}[ds} Sj (t')s k (t")(eMW))f , (43) 

and similarly 

S Sj (t')/Se k (t") = ( Sj (t')s k (t")) f = JI [ds}[ds} Sj (t')s k (t''){eMW))f . (44) 

It is not obvious how to calculate (exp(W))f directly, since W involves terms like JijSiSj 
which couple different sites. Following Sompolinsky and Zippelius [76], and Narayan and 
Fisher [19], we can circumvent this problem by performing a change of variables from the 
spins Sj to local fields Jfji = Y.j JijSj- (We introduce the coefficient J on the left hand side to 
keep the dimensions right.) At the saddle point of the associated action the new variables fjj 
(for all j) are given by the mean- field magnetization and the different sites become decoupled. 
A saddle point expansion becomes possible, because the coefficients in the expansion can be 
calculated in mean- field theory — they are also the same for all sites j. 

Here is how it works: we insert into Z the expression 

1 = 1/2tt J J[df}][df}]J[rj\exp(iY, J " £ J ^ J Ut))) (45) 

j i 

where J^[fj] stands for the suitable Jacobian, which is simply a constant and will be sup- 
pressed henceforth. Integrating out the auxiliary fields fjj, one recovers that the expres- 
sion in eq. (45) is the integral over a product of (^-functions which impose the definitions 
Jfjiit) = J2j Jij s jif) at all times for all i. 

After some reshuffling of terms and introducing some redefinitions that are motivated by 
the attempt to separate the nonlocal from the local terms, one obtains 

Z= [ [ [d^mUZjifjj^eM- I dtfj^C^J^Jm)} (46) 
j J I 
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where Zj[fjj,fj 3 ] is a local functional 



z AVj,Vj}=J J[ds j ][ds j ]{ex.p{J 1 J dt[Jrj j (t)s j (t)+ 
is 3 (t)(d t s 3 (t)/T - Jf, 3 {t) -H-f j + 5V/5 Sj )]})f 



(47) 



(we have absorbed a factor % in the definition of 77 ). In short this can also be written as 



with the effective action S e ff, now expressed in terms of the "local field" variables fj and r/ 



as the random field averaged probability distribution for the possible paths the system can 
take through the spin configuration space as the external magnetic field is slowly increased. 
Each path is specified by a set of N effective field functions Jrji(t) (i runs over the lattice 
with iV — > 00 spins). Z is the integral of this (normalized) probability distribution over all 
possible paths of the system and is therefore equal to I. 

The stationary point [fj^Vj] °f the effective action is given by 




(48) 




(49) 



Physically we can interpret the functional 




(50) 



[ SS eff/Srjj\ip fj° = ' 
''3' h 



(51) 



and 



[SSeff/Sfjjj^o = . 



(52) 



With eqs. (47) and (49) we find the saddle-point equations: 



(53) 



3 



and 
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- E ^? = • (54) 

3 

The notation () -o „ here denotes a local average, obtained from the local partition function 
Zi, after having fixed f\ i and fji to their stationary-point solutions ff* and ff^. 
For example (from eqs. (47), (49), (52), and (54)) 

f / j j i 

(expjj- 1 1 dtE^W^W 

i 

+i^(0(9 tSj (0/r - Jig(f) -H-f j + 5V/5 Sj )]}) f . (55) 
Eq. (53) and eq. (54) have the self-consistent solution 

%(t) = 0, (56) 

and 

= M(t) = , (57) 

where M(t) is the random field average of the solution of the mean field equation of motion 

d tSj (t)/r = Jffl(t) +H + f,- SV/Ssj . (58) 

(This can be seen by setting rj° = in eq. (55). Integrating out the s fields we see that 
Zj is the random field average over a product of ^-functions which impose eq. (58) by their 
argument. Eq. (57) is the self-consistency condition for this mean-field equation of motion.) 

We can now expand the effective action S e g in the variables fjj = (fjj — 77°) and rjj = 
(fjj — fjj), which corresponds to an expansion around mean-field-theory: 

Z = f J[d V ][dri]eMSeff) (59) 
with an effective action (expressed in the new variables rj and fj): 



dt 



m+n 



r 1 /' 

Se ff = -E du^j^tMt) +EEznj* 

j,l j m,n=0 ■ ■ 

(t m )r)j(t m+1 ) ■ ■ ■ r)j(t m+n ) (60) 
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Here, as seen by inspection from eq. (49) and eq. (47) 



d 



d 



U r 



df)j(t m+1 ) dr]j(t m+n ) 
d d 



Sfjjih) ■ ■ ■ 5f)j(t 



(61) 



r]=0,r)=0 



{(8(h) - Tl^h)) ■ ■ ■ (S(t m ) - V°(.t m )))l,rp,rp , 



de(t m+1 ) de(t m+n ) 

i.e., the coefficients u mn are equal to the local (1), connected responses and correlations in 
mean field theory. 

Again, local (/) means [19] that we do not vary the local field rjj in the mean-field equation 



-Ua,-(f) = Jrg(t) +H + fj- ^ + Je{t) 



(62) 



when we perturb with the infinitesimal force Je(t). 



B. Source terms 

Correlations of s and s can be related to correlations of r\ and fj [19]. If we introduce the 
source terms 

J dtisjMijW-isjitMt)) (63) 

into the action, we can write the correlations of s and s as functional derivatives with respect 
to I and I at I — I — 0. A shift in the variables rj and 17 by I and I respectively leads to a 
source term of the kind 

jj-^m-kmvAv-iAt)) (64) 

so that derivatives with respect to / and 1 give correlation functions of 57 and 77. For low 
momentum behavior the factor JJ^ 1 can be replaced by one since J2i Jij 1 = J 1 - 



C. Implementing the history 



Up to here the effective action S e g manifestly involves the entire magnetic field range 
—00 < H < +00. As we discuss in appendix C it turns out, however, that in the adiabatic 
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limit a separation of time scales emerges. The relaxation rate kT in response to a pertur- 
bation is fast compared to the driving rate Q/k of the external magnetic field. The static 
critical exponents can then be extracted self-consistently from a RG calculation performed 
at a single, fixed value H of the external magnetic field. The analysis is much simpler 
than one might have expected. Instead of dealing with the entire effective action which 
involves all field values H, it suffices in the adiabatic limit to calculate all coefficients u mn in 
eq. (60) at one single fixed magnetic field H, and then to coarse grain the resulting action 
S e ff(H) = Sh- There are no corrections from earlier values of the external magnetic field. 

Physically this corresponds to the statement that increasing the magnetic field within 
an infinite ranged model (mean field theory) and then tuning the elastic coupling to a short 
ranged form (RG) would be equivalent to the physical relevant critical behavior, which 
actually corresponds to first tuning the elastic coupling to a short ranged form and then 
increasing the force within a short ranged model [20] . In their related calculation for CDWs 
below the depinning threshold [20], Narayan and Middleton give an argument that this 
approach is self-consistent for their problem. In the appendix C we first show that their 
argument applies to our system as well, and then discuss the consistency of the magnetic 
field decoupling within the RG treatment of the entire history for separated time scales. 

Where did the history dependence go? Note that the values of the coefficients u mn at 
field H are still history dependent (in the way the mean-field solution is). Also, causality 
must be observed by the coarse graining transformation, so that even in the adiabatic limit 
the intrinsic history dependence of the problem does not get lost. 

D. Calculating some of the u mn coefficients at field H 

In appendix C we show that u mn basically assume their static values in the adiabatic 
limit. In this section we will briefly outline their derivation and quote the relevant results. 

We have to be consistent with the history of an increasing external magnetic field, when 
expanding around the "mean-field-path" rf(t). This implies that for calculating responses 
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from eq. (61) we must only allow a perturbing force Je(t) that increases with time in eq. (62). 
For example, for u 1:1 we add a force Je(t) = JeQ(t — t') in eq. (62), where Q(t — t') is the 
step function, and solve for (s(t)\H+j e (t))f- The local response function is then given by 
the derivative of lim t ^o[(((s(t)\ H +je(t)} / — ( s (^)|//)/)/ e ] with respect to (—£'). The higher 
response functions are calculated correspondingly. (The most important ones are calculated 
in appendix C.) One obtains in the low frequency approximation for the first few terms of 
the effective action of eq. (60) at field H 



with 



and 



Sh = ~J:J dttfjTj&Mt) -£ I dt^(t)[-ad t /T - uff] Vj (t) 

3,1 3 

+ W di>fe(t)Mt)) 3 

j J 

+ J2j dtl J dt2 | M 2,0%(^l)%(^2) , (65) 



a = (j/k + Api-Jr] - H + k))/k, (66) 

uff = 2Jp(-Jr]° - H + k) + J/k, (67) 

w = -2J 2 p'{-Jrf -H + k), (68) 

u = 2J 3 p"(-Jr]° -H + k), (69) 



f -H-r]°+k r-H-rp+k 
U 2fl = R 2 /k 2 + 4( / p{h)dh)-A{ p(h)dh) 2 

p-H-rf+k 

-4( / (h/k)p(h)dh) . (70) 

Eq. (70) implies that -u 2 ,o > for any normalized distribution p(f). 
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VII. PERTURB ATIVE EXPANSION 



A. The Gaussian theory for d > d c : response and correlation functions 

One can show [53,66] that for systems with dimension d above the upper critical dimen- 
sion d c , the nonquadratic terms in the action become less and less important on longer and 
longer length (and time) scales. Near the critical point, where the behavior is dominated 
by fluctuations on long length scales, the system is then well described by the quadratic 
parts of the action, and the calculation of correlation and response functions amounts to 
the relatively simple task of solving Gaussian integrals. It should come as no surprise that 
the mean-field exponents are recovered, since the quadratic parts of the action represent the 
lowest order terms in the saddle-point expansion around mean-field theory. 

In our problem the action Sh of eq. (65) has the quadratic part 



Q(v,v) = -£/ dtJ^jfijitMt)-^ J dtf) j (t)[-<ft/r -u$*\v s (t) 

3,1 3 

+(1/2) I dh [ dt 2 %(*i)%(*2Ko . (71) 

3 J J 

In the long-wavelength limit we can write J~ l (q) — 1/J+ J^q 2 [64]. 
Rescaling fj, u and q we can replace the constants J2J and a by 1. The low frequency 
part of the 7777-term in Q(i], 17) is then given by 

- Jd d qJ dtrji-qM-dt/ro + q'-X-'/JHQit) (72) 

where 

x - 1 = J( u ff - 1) = 2J 2 p{-JM - H + k) - J{k - J)/k (73) 

is the negative static response to a monotonically increasing external magnetic field, calcu- 
lated in mean-field-theory. 15 



15 Note that % 1 = J(k — J)/kt, where t is a parameter used in the mean-field scaling functions 
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In frequency space, it can be written as 

- J d d g J dt(-iu/r + q 2 - x-'f-m-q, -w)v(q, ■ (74) 

The fjfj term in Q(rj,fj), given by 

J d d q J dh J dt 2 l/2u 2fi f}{-qM)fl{qM) (75) 
can also be expressed in frequency space 

J d d q J du(l/2)u2, 8(u)r)(-q, -u)r)(q, u) . (76) 
The expressions from eq. (74) and eq. (76) can be written together as 

Q(V,V) = - J dw J d d q (77) 

(rj(-q,-u>), ni-q-uj)) ( -l/2u 2 ,o5(w) (iui /F +q 2 - X _1 / 'J) \ ( vil,") \ 

\ (-iu/ToW-X-i/J) J \ V (q,u) ) 

This can be used to determine the response and correlation functions at field H at low 
frequencies to lowest order in the expansion around mean field theory [77]. Inverting the 
matrix one obtains [19] 

Gr, v (q, u) = {r)(-q, -u)r}(q, u)) w l/(-iw/r + q 2 - X^/J) (78) 

and 

G m (q, u) = (v(~q, -u)ri(q, u)) w u 2fl 5(u))/\ - iuj/T + q 2 - x" 1 / J\ 2 ■ (79) 

The 5{uj) in eq. (79) is a consequence of the underlying separation of time scales. It will 
lead to an essentially static character of the RG analysis of the problem. This might have 
been expected, since the critical phenomena we set out to describe are essentially static 



in appendix A. It is defined for the soft spin model in eq. (A29). In particular, this implies that 
near the critical point \~ l scales with the original parameters r and h in the same way as t, see 
eq. (A 4) in the same appendix. 
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in nature. At the critical point x^ 1 = we have Gjj V (q,cu = 0) ~ q 2 ~ v with r\ = and 
G m (q,u = 0) ~ g 4- *' with fj = at lowest order in perturbation theory. One can Fourier 
transform the correlation functions back to time 



and 



Gfj V (q,t,t') = (rj(q,t)r](-q,t')) 

r exp{-r (g 2 - X^/JW - t)} for t' > t 
for t' < t 



= < 



G m (q, t, t') = J dti J dt 2 Gf jT] (q, ~(t' -t)+t 1 + t 2 )u 2 ,oGfj V (q, t 2 ) 
= u 2 , /(q 2 -X- 1 /J) 2 - 



(80) 



(81) 



B. The RG analysis 

In dimension d < d c the nonquadratic parts of the action are no longer negligible near 
the critical point. One obtains corrections to the mean-field behavior. 

The Wilson-Fisher coarse graining procedure is an iterative transformation to calculate 
the effective action for the long wavelength and low frequency degrees of freedom of the 
system. In each coarse graining step [78,65,19,66] one integrates out high momentum modes 
of all frequencies f)(q,u) and rj(q,uj), with q in a momentum shell [A/6, A], b > 1, and 
afterwards rescales according to q = fe~V) u) = 5 _1 u/ ', r)(q,uj) = b Cp fj'(q,uj), and r](q,u) = 
b Cp rj'(q,u). As usual, the field rescalings c p and c p are chosen such that the quadratic parts 
of the action at the critical point (x -1 = 0) remain unchanged, so that the rescaling of the 
response and the cluster correlation function under coarse graining immediately gives their 
respective power law dependence on momentum (i.e. this is an appropriate choice of the 
scaling units.) Without loop corrections this implies that z = 2, r)(x,t) = b~ d / 2 ~ z fj'(x,t) and 

Tj(x,t) =b- d ' 2+2 i(x,t). 

Performing one coarse graining step for the expansion for Sh or eq. (65) yields a coarse 
grained action which can be written in the original form, with "renormalized" vertices u' mn . 
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Without loop corrections, the vertices of the coarse-grained action are simply rescalings 
of the original vertices, which can be easily read off using the rescalings of q, u, r], and 
r). Taking into account that each 5/5e(t) involves a derivative with respect to time, and 
therefore another factor b~ z under rescaling, we arrive at 

<„ = 6 [ - (m+n)+2ld/2+2 ^ m ,„ • (82) 

This shows that above 8 dimensions all vertices that are coefficients of terms of higher 
than quadratic order in the fields, shrink to zero under coarse-graining and are therefore 
"irrelevant" for the critical behavior on long length scales and at low frequencies. 

The "mass" term x" 1 in the action actually grows under rescaling in any dimension: 

= b'ix- 1 ) (83) 

(without loop corrections). At — an d if all irrelevant coefficients are set to zero, 
the action does not change under coarse graining. This case is obviously a fixed point of 
the coarse graining transformation, where the systems looks the same on all length scales. 
There is no finite (correlation) length determining the long wavelength behavior, which is 
just what one would expect at a critical point. 

The critical exponents can be extracted from a linearization of the transformation around 
the associated fixed point. In the RG sense "relevant" eigenvectors of the linearized transfor- 
mation (i.e. coefficients like x _1 that grow under coarse graining) render the corresponding 
scaling fields or tunable parameters. 

Below eight dimensions the vertex -u 1)2 = w — 2Jp'(—Jrf — H + k) is the first coefficient 
of a nonquadratic term to become relevant. An action with the original parameters x 1 = 
and w 7^ corresponds to a system with less than critical randomness R < R c at the 
onset field of the infinite avalanche. In appendix E we show how to extract the mean- 
field exponents for the infinite avalanche line from the scaling above eight dimensions and 
that the RG treatment suggests a first order transition for the same systems in less than 8 
dimensions. 
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In systems where the bare value of w is zero at the critical fixed point with x~ l — 0, 
all nonquadratic terms are irrelevant above six dimensions. As can be seen from appendix 
A, this case constitutes the interesting "critical endpoint" at R = R c and H = H C (R C ), 
which we have discussed in the introduction. The corresponding "upper critical dimension" 
is d c = 6. For any dimension higher than 6, the mean-field exponents should describe 
the behavior on long length scales correctly [53]. For d < 6, the vertex = u becomes 
relevant, while all higher vertices remain irrelevant. We are left with the effective action 
which includes all vertices relevant for an expansion around 6 dimensions: 

S=- jd d qjdt fj(-q, t) [-dt/Vo + q 2 - x" 1 / J] r)(q, t) (84) 
+(1/6) E / dtf b (t)(r lj (t)fu+ (1/2)^ J dh J dt 2 t^Ot^Ko 

3 3 

Below 6 dimensions, where u does not scale to zero, we will perform the coarse graining 
transformation in perturbation theory in u. At the fixed point, in 6 — e dimensions, u will 
be of 0(e). The perturbation series for the paramenters in the action and thus also for 
the critical exponents, becomes a perturbation series in powers of e. From the form of the 
action one can derive Feynman rules (see appendix C), which enable us to write down the 
perturbative corrections in a systematic scheme. Examples of their derivation for the 4 
model are given elsewhere [66,53,64]. 

1. Loop corrections 

In the remaining parts of this section we perform a coarse graining transformation to 
first order in e. 16 From the integration over the short wavelength degrees of freedom (of all 



16 In principle this is not necessary, since we can read off the corrections to 0(e 5 ) from the mapping 
to the e-expansion of the pure Ising model in two lower dimensions, which we discuss in the next 
section. However we will need the techniques introduced here later for the calculation of the 
avalanche exponents. 
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frequencies) one obtains loop corrections to various vertices. Figure 11 (a) and (b) shows the 
corrections to x 1 an d u which are important for an 0(e) calculation. The dots correspond 
to the vertices with the indicated names. An outgoing arrow corresponds to an fj operator, 
an incoming arrow corresponds to an r\ operator. 

We consider the fjrj term in the action as propagator and all other terms as vertices. An 
internal line in a diagram corresponds to the contraction 



r exp{-r (g 2 - x-'/JW - t)} iort>>t 

(85) 

for t' < t 



with q in the infinitesimal momentum shell A/b < q < A (b > 1) over which is integrated. 
This expression can be approximated by 5(t — t') in the low frequency approximation [19]. 
Note, however, that causality must be obeyed, i.e. t' > t. Figure 11 (c) shows an example 
of a diagram that violates causality and is therefore forbidden. External (loose) ends in 
a diagram correspond to operators that are not integrated out, i.e. modes of momentum 
q < A/b outside of the momentum shell. Each internal line carrying momentum contributes 
a factor 

l/(g 2 - X - l /J) ■ (86) 
The entire loop in diagram 11 (a) contributes the integral 

h= £ d d q/(2nri/(q"-x' 1 /J) 2 (87) 
J A/b 



(integration over time is already performed). Similarly the loop diagram in figure 11 (b) 
yields the integral 



h= Cd d q/(2 7 r) d l/(q 2 - X - 1 /J) 3 . (88) 

J A/b 



After each integration step we also have to rescale momenta, frequencies and fields. Thus 
the recursion relations for x 1 1 J an d u become 

(x~ 1 /J)' = b 2 {x' 1 /J+^Qh) (89) 
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and 



u' = b e - + 



3! 2! 



2 



36/ 2 . (90) 



(There are no loop corrections to U2,o at this order, so u' 20 = W2,o-) The integrals I\ and J 2 
have to be computed in 6 — e dimensions. To that end one uses the relation [53] 

I = J (&/(^) = WW I dq^M) , (91) 

where Sd is the surface area of a unit sphere in d dimensions: 

S d = 27r d/2 /T(d/2) . (92) 

In performing this analytic continuation in dimension we make the very strong assumption 
that the physical properties vary smoothly with dimension. It has been justified by its great 
success in many problems (although it gives an asymptotic expansion in e which is not well 
behaved, see appendix D). Disregarding potential complications we then expand [53] 

/i, 2 = /i, 2 (0) + < 2 (0) + O(e 2 ). (93) 

Since both I\ and J 2 are multiplied by u in the recursion relations, which will be of 0(e) at 
the fixed point, we only need -^1,2(0) for an 0(e) calculation. With Kd = Sd/ (2n) d we obtain 
[53] 



h = [ A d d q/(2n) d l/(q 2 - X ~ l / jf 
Jk/b 

= d d q/(2n) d (l/ q ')(l + 2( X ' 1 /J)/q 2 + O^ 1 / j f)) 

= K 6 A 2 (1 - l/b 2 )/2 + 2K 6 ( X - 1 /J) In 6 + O^ 1 / J) 2 , e) . (94) 



Similarly 



I 2 = K Q \nb + 0{{ X - l /J)\t)- (95) 
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2. Recursion relations to O(e) 



With v = M2,oW and b e = e elnb = 1 + elnfc + 0(e 2 ) the recursion relations then become 
(^7 J)' = b\ X - l /J + *V(47r) 3 A 2 (l - l/6 2 )/4 + v/^fix^/J) In ft) (96) 



and 



v' = v + v[e + 3w/(47r) 3 t;] Info. 



(97) 



There are two fixed points of these relations. The Gaussian fixed point (% V J)* = 0,v* = 0, 
and a new, nontrivial fixed point at 



(x-'/jy = -eA 2 /12 



(98) 



and 



v * = -(47r) 3 e/3, 



(99) 



which is sometimes called the Wilson-Fisher (WF) fixed point. We linearize the recursion 
relations around these fixed points to calculate the critical exponents. At the Gaussian fixed 
point the corresponding matrix is 



V 



/ 



V 



o 



d {*r) Id bC 1 1 J) dix-'/jy/dv 
dv'/d(x~ 1 /J) dv'/dv 
At the Wilson-Fisher fixed point it is 

V-e/3) 3A!(fo2_l)/( 47r )3\ 



/ 6 2 f (6 2 - l)/(4vr) 3 X 



(100) 



/ 



V 



(101) 



b' e 

The eigenvalues of the transformation linearized around the Gaussian fixed point are 
given by e x -\ = b yt with y t = 2 and e 2 = 6 e with eigendirections e x -\ = (1, 0) and 
e 2 = (A 2 /(4(4vr) 3 ), 1) in the (x -1 /J>) plane. 

The eigendirections of the transformation linearized around the Wilson-Fisher fixed point 
are e x -\ = (1,0) and e 2 = (— eA 2 /12, — (47r) 3 e/3). The corresponding eigenvalues, called 
A( x -i) and A 2 are given by 
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(102) 



with y t = 2 — e/3 and 



A? = 6" 



(103) 



5. Flows and exponents 

a. For e < 0, i.e. d > 6: the Gaussian fixed point is unstable along the e x -i direction, 
but stable along the e 2 direction. Thus v (or equivalently u) is an irrelevant variable. 
The correlation length scales as Hx^/J) — ^((x -1 /^)') since under coarse graining the 
coordinates scale as x' = bx. From the appendices A and C we know that w = implies 
that h = and ~ r. 17 One obtains (along h = 0) 

f (r) ~ ^(fi^r) . (104) 

Choosing 6 = r" 1 / 2 " it follows that 

£ = r -iM£ (1) . (105 ) 
Since ^ ~ r _I/ along h = 0, one obtains 

i/ = l/y t . (106) 

Therefore v = 1/2 for d > 6. 

6. For e > 0, i.e. d < 6: the Gaussian fixed point becomes unstable along both e x -\ 
and e 2 . Along the e 2 direction the flow is towards the attractive WF fixed point, which 
will govern the critical behavior of all systems that do not initially sit at the Gaussian fixed 
point u = Jp"(—JM — H) = 0. For those systems one then obtains (see eq. (102)) 



17 See appendix A, eq. (A13) in particular, keeping in mind that the parameter t given there is 
proportional to x" 1 -, apart from some irrelevant adjustments for the soft-spin model. 
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v = l/y t = 1/2 + e/12 (107) 

The corresponding flows in the (x~ 1 /J^ v ) plane are shown in figure 12. 

As in the RG treatment of the standard Ising model, the anomalous dimension rj is equal 
to zero to first order in e, since the field rescalings acquire no loop corrections to first order 
in e [53]. (The coefficient of the g 2 -term in the propagator and the i^o-term only receive 
0(e 2 ) and higher order loop corrections.) Since the iu term also receives no corrections to 
0(e), we have z = 2 + 0(e 2 ). 

4- Equation of state 

We can also calculate 0(e) corrections to the equation of state [79], and consequently for 
the entire scaling function of the magnetization in eq. (9). Following Wallace [79] we replace 
the field fj(x) in the Hamiltonian by a new field L(x) with zero expectation value 

f)(x)=L(x)+M (108) 

with 

(L) = 0. (109) 

In other words, in contrast to the preceding analysis, here we expand the action around 
the true, not mean-field, magnetization M at the considered order. At tree level [L(x) = 
r)( x )]tree-ievei and [M = r)°]tree-ievei and [u lfi (M) = % ree -ievei- At higher order however, unlike 
in an expansion around mean-field theory, ui^(M) is no longer zero. Its role is to cancel 
off all tadpole graphs [79,19], so that (L) = at any given order. The lowest order tadpole 
correction is shown in figure 13. Here we follow Wallace's calculation in [79]. Instead of 
performing an iterative integration over successive momentum shells, as we have done before, 
all momentum modes are taken into account in one step in which the integration extends 
from up to the cutoff A [68]. The analytic expression to be added to w^o, corresponding 
to the lowest order correction of figure 13, is then given by 
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(l/2)u 2fi w [ l/(q 2 - X - l ) 2 d d q. (110) 
As we have explained, the true magnetization M is chosen such that 

«i,o(M) + (l/2)u 2 , w [ A l/(q 2 - x^f^q + • • • = (111) 

J 

where • • • stands for all higher order corrections. Eq. (Ill) then constitutes the equation 
of state at the given order. In an expansion around rf = M (where M is now the true 
magnetization), u lj0 (M) is given by eq. (61), eq. (62) and section VIC to be 

r -JM-H+k 

u 1>0 (M) = (k + H)/k-(k-J)M/k-2 p(h)dh. (112) 

J —oo 

Close to the critical point the integral can be expanded around —JM—H+k = —Jm — h = 0, 
where m = M — M c is the deviation of the magnetization M from the mean-field critical 
value M c = +1, and h = H — H C (R C ), with the mean- field value H C (R C ) = k — J, and 
r = (R c - R)/R with R c = 2kJ/(V2ir(k - J)). One obtains 

«i,o = \ + y™ + ^P"(0)m 3 + • • • , (113) 

where p"(h) denotes the second derivative of the distribution of random fields with respect to 
its argument and • • • stands for higher order corrections. To calculate the equation of state at 

one loop order we insert into eq. (Ill) the expressions x^ 1 — —2Jr/ (y/2nR) + \/2um 2 H 

(see eq. (73)), w = urn + ■ ■ ■ (see eq. (68)) and u = 2J 3 p"(0) + • • • (see eq. (69)), where • • • 
denotes higher orders in m and h. One then finds that the calculation is analogous to the 
one done for the Ising model in two lower dimensions. For details on solving the integral 
etc. we refer the reader to the article on the equation of state for the Ising model in 4 — e 
dimensions by D.J. Wallace in reference [79], especially equation (3.35). (In fact, with the 
following formal identifications, the resulting equations of state in the two systems can be 
mapped onto each other: h/k = —h w , 2Jr / '(y/2nR) = —t w , 2J 3 p"(0) = —(uo) w , m = —m w . 
We have denoted the quantities in Wallace's article by an index "w".) 
We find for the equation of state to one loop order: 
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h I 9 7 3 \ , 

Ul o = = - + 2J/(V2nR)rm + — p"(0)m 3 + -em{2J/(V2^R)r 
k 3! 6 

+J 3 p"(0)m 2 ) \n(-2J/(V^R)r - J 3 p" '(0)m 2 ) . (114) 
Here e = 6 — d and r and m are the reduced randomness and magnetization respectively. 

5. The exponents (3 and 5 

This equation allows us to extract the critical exponents f3 and 5 by taking the cor- 
responding limits: The exponent 5 can be obtained by setting r = 0. That leaves one 
with 

= ui >0 = h/k + l/3J 3 p"(0)m 3+e . (115) 
Since by definition of 5 one has \h\ ~ \m\ 5 , we find 

5 = 3 + e. (116) 
Similarly one sets h — to get /? and finds from |r| ~ Imj 1 /' 3 (for i? < R c ) that 

"-5-5- (117) 



6 1 . The entire scaling function 

From eqs. (114), (116) and (117) one can construct the entire scaling function /(r/m 1///3 ) 
to 0(e) in 

\h/k\ = \m\ s f(r/\m\ 1/p ) (118) 
directly, as shown in reference [79]. With the redefinitions 

m 2 new = -(l/3)J 3 p"(0)m 2 (119) 

and 
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x = -(2Jr)/(v / 2Ti?|m new | 1 /' 3 ) (120) 

(for r > 0), one obtains 

f(x) = x + 1 + l((x + 3) ln(x + 3) - 3(x + 1) ln(3) + 2x ln(2)) + 0(e 2 )) , (121) 

in 

\h/k\ = \m new \ 5 f(x). (122) 

In the next section we discuss a formal mapping of our e expansion to the e expansion of 
the pure Ising model in two lower dimensions, which allows us to copy the 0(e 2 ) correction 
to the equation of state from previous calculations for the pure Ising model [80]. 

VIII. MAPPING TO THE THERMAL RANDOM FIELD ISING MODEL 
A. Perturbative mapping and dimensional reduction 

We notice that the results found for u, (5 and 5 in 6 — e dimensions (eqs. (107), (117), and 
(116)) are the same as those for the regular equilibrium Ising model in 4 — e dimensions. In 
fact, this equivalence actually extends to all orders in e: We will show that the e-expansion 
for our model is the same as the e-expansion for the equilibrium random field Ising model 
to all orders in e [130]. Once this equivalence is established, we can use that the (6 — e)- 
expansion of the equilibrium random field Ising has been mapped to all orders in e to the 
corresponding expansion of the regular Ising model in two lower dimensions [81,82]. 

The easiest way to recognize that the e-expansion for our model and for the equilibrium 
RFIM should really be the same is by comparing the corresponding effective actions. In a 
dynamical description of the equilibrium RFIM at zero external magnetic field the following 
effective Langevin equation of motion for the spin-field 0(r, t) was used [60] 

d t (f)(x, t) = -r o (-V 2 0(:r, f) + r 0( j)(x, t) + 1/Qg (f) 3 (x, t) - h R (x)) . (123) 
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Iir(x) represents spatially uncorrelated quenched random fields distributed according to a 
Gaussian of width A and mean zero. hx(r, t) is the thermal noise field, which is taken to 
be gaussian with vanishing mean value and the variance 

(h T (r, t)h T (r', t')) = 2kT/T 5(r - r')5(t - t') . (f 24) 

The corresponding Martin Siggia Rose generating functional is 

thermal = J ^ J [d(f)] exp{ _ J ^ J ^ ( _^ t )(-d t /T + q 2 + r )<f>(q, t)) 

+ Jd d xJ dt(4>(x,t)(-l/6g )(f>(x,t) 3 ) 

+ Jd d xJ dt(i(x,t)(h R (x) + hr(r,t))}. (125) 

Since again 2j^ rma ^ — 1, we can average the partition function directly over the random 
fields h R and the thermal noise h T : 

rythermal I ry thermal \ ( ~[ oa\ 

Z H - { Z H )h R - U2bJ 

The average over the random fields at each x and over the thermal noise fields at each x 
and t yields, (after completing the square): 

2 thermal = J J {- j d « q j dt 0(-g, t)(-d t /T + q 2 + r )<j>(q, t) 

+ J d d x J dt 4>{x,t){-\g G )<t?{x,t) 

+ J d d x j dti J dt 2 (j)(x,t 1 )(j)(x,t 2 )A 2 /2 

+ J d d x J dt(j) 2 (x,t)(2kT)/T ] j . (127) 

With the identifications r = — x _1 (iJ ), u = —l/6g , u 2 ,o = A 2 and T = 0, we see 
that the argument of the exponential function is the same action as the effective action for 
our zero temperature, nonequilibrium model in eq. (65). Setting T to zero in the action for 
the equilibrium RFIM does not change the expansion for the static behavior, since it turns 
out that corrections involving temperature are negligible compared to those involving the 
random magnetic field [82,60,83] — the temperature dependence is irrelevant in the thermal 
RFIM and the time dependence is irrelevant in our zero-temperature dynamical RFIM, 

52 



leaving us with the same starting point for the calculation! This equivalence implies that 
in 6 — e dimensions we should obtain the same critical exponents to all orders in e for our 
model, as were calculated for the thermal random field Ising model, which in turn are the 
same as those of the pure equilibrium Ising model in 4 — e dimensions [81,82]. 

This observation is rather convenient, since it provides us with results from the regular 
Ising model to 0(e 5 ) for free. In 6 — e dimensions we read off [59] 

\jv = 2 - e/3 - 0.1173e 2 + 0.1245e 3 - 0.307e 4 + 0.951e 5 + 0(e 6 ) (128) 

7] = 0.0185185e 2 + 0.01869e 3 - 0.00832876e 4 + 0.02566e 5 + 0(e 6 ) (129) 
/? = 1/2 — e/6 + 0.00617685e 2 - 0.035198e 3 + 0.0795387e 4 

-0.246111e 5 + O(e 6 ) (130) 
(35 = 3/2 + 0.0833454e 2 - 0.0841566e 3 + 0.223194e 4 - 0.69259e 5 

+0(e 6 ). (131) 

(3 and 5 have been calculated from r\ and v using the perturbative relations: (3 = (u/2)(d — 
4 + fj) and 5 = (d — 2rj + fj) j(d — 4 + fj) [84] , with rj = fj to all orders in e [60] . 

By the same mapping we obtain the universal scaling function for the magnetization to 
0(e 2 ) [80,79]. We set 

h = m 5 f (x = r/m 1/p ) (132) 
in which the normalizations of x and the universal scaling function / are chosen such that 

/(0) = 1, /(-1) = 0. (133) 
The expansion to second order in e is then 

f{x) = l + x + ef 1 (x) + e 2 f 2 (x) (134) 

with 



= l[(x + 3) ln(x + 3) - 3(x + 1) ln(3) + 2x ln(2)] 
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(135) 



and 



18 



{[6 In 2 -91n3][3(x + 3) ln(x + 3) + 6x\n2 - 9(x + 1) In 3] 



+ §(:r + l)[ln 2 (:r + 3)-ln 2 3] 

+36 [In 2 (x + 3) - (:r + 1) In 2 3 + adn 2 2] 



541n2[ln(x + 3) + x\n2- (x + l)ln3] 



+25[(x + 3) ln(x + 3) + 2xln2 - 3(x + 1) In 3]} . 



(136) 



The scaling function f(x) has actually been calculated up to order e 3 [67]. As it stands the 
expression (134) meets the Griffith analyticity requirements [85] only within the framework 
of the e-expansion, but not explicitly. These subtleties can be avoided by writing it in a 
parametric form [80], which can then be compared directly with our numerical results for 
the universal scaling function of dM/dH in 5, 4, 3, and 2 dimensions. We will present the 
results in a forthcoming paper [45]. 

The dynamic exponent z cannot be extracted from the mapping to the regular Ising 
model. It was calculated separately to 0(e 3 ) for the equilibrium RFIM [60] and found to be 
given to this order by 1819 



18 Eq. (137) is only a perturbative result for z which does not reveal the presence of diverging barrier 
heights that lead to the observed slow relaxation towards equilibrium [2,62,3-5]. Nonperturbative 
corrections are expected to be important in the equilibrium random field Ising model. 

19 While we expect the 6 — e results for the static exponents 0, 5, u, rj, fj, r, a, etc. to agree with 
our hard-spin simulation results close to 6 dimensions, this is not necessarily so for the dynamical 
exponent z. There are precedences for the dynamics being sensitive to the exact shape of the 
potential, at least in mean-field theory. Such differences have been found for example in the case 
of charge density waves [17,19,20,127], if studied for smooth potentials, for cusp-potentials, or for 
sawtooth potentials. 



z = 2 + 2r] = 2 + 0.037037e 2 + 0.03738e 3 + 0(e 4 ) . 



(137) 
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Because of the perturbative mapping of our model to the equilibrium RFIM, eq. (137) also 
gives the result for z in our nonequilibrium hysteretic system. 

We have performed a Borel resummation (see also appendix D) of the corrections to 
0(e 5 ) for 77, v, and the derived [84] exponent (35. The exponent (3 is then given through 
(3 = (35 — (2 — r\)v (eq. (24)). Figure 14 shows a comparison with our numerical results in 3, 
4, and 5 dimensions. 

The agreement is rather good near 6 dimensions. However the apparent dimensional 
reduction through the perturbative mapping to the Ising exponents in 2 lower dimensions 
gradually loses its validity at lower and lower dimensions. It is after all only due to the 
equivalence of two asymptotic series, both of which have radius of convergence zero. Table I 
shows a comparison between the numerical exponents for our model and for the equilibrium 
RFIM in three dimensions. 

B. Nonperturbative corrections 

The mapping of the e-expansion for the thermal random field Ising model to the expan- 
sion for the Ising model in two lower dimensions has caused much controversy when first 
discovered. The problem was that it had to break down at the lower critical dimension, 
where the transition disappears. There is no transition in the pure Ising model in d = 1, 
but the equilibrium RFIM is known rigorously to have a transition in d = 3 [102,103]. The 
same is true for our model: numerical simulations indicate [45,98] that the lower critical 
dimension is lower than three — probably equal to two. 

In the case of the equilibrium random field Ising model it was finally agreed [62,82] that 
this breakdown might be due to nonperturbative corrections. The point is that proving a 
relation to all orders in e does not make it true. In the equilibrium RFIM there are at least 
two sources of nonperturbative corrections: 

a. The "embarrassing" correction: It was found that there was a calculational error 
in the (6 — e)-expansion for the RFIM. The perturbation series was tracing over many 
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unphysical metastable states of the system, instead of just taking into account the ground 
state, which the system occupies in equilibrium. There were indications that this error leads 
to nonperturbative corrections, which would destroy the dimensional reduction outside of 
perturbation theory [82]. 

In our calculation we have avoided the embarrassing source of nonperturbative corrections 
found in the equilibrium random field Ising problem. Given the initial conditions and a 
history H(t), the set of coupled equations of motion for all spins will have only one solution. 
In the Martin Siggia Rose formalism, the physical state is selected as the only solution that 
obeys causality, there are no unphysical metastable states coming in. Therefore we believe 
our results should also apply to systems below the critical randomness, at least before the 
onset of the infinite avalanche. 

b. Instanton corrections: Even without the embarrassing correction, there is no reason 
why a perturbative mapping of the expansions about the upper critical dimensions should 
lead to a mapping of the lower critical dimensions also. The e-expansion is only an asymp- 
totic expansion — it has radius of convergence zero. As we discuss in appendix D, there is 
no known reason to assume that the e-expansion uniquely determines an underlying func- 
tion. It leaves room for functions subdominant to the asymptotic power series: If the series 
So° fk zk is asymptotic to some function f(z) in the complex plane as z — > 0, then it is also 
asymptotic to any function which differs from f(z) by a function g{z) that tends to zero 
more rapidly than all powers of z as z — > [90] . An example of such a subdominant function 
would be g{z) = exp(— 1/z). While some asymptotic expansions can be proven to uniquely 
define the underlying function, this has not been shown for the e-expansion (appendix D) 
- not for our problem, nor for the equilibrium pure Ising model, nor for the equilibrium 
thermal random field Ising model. 

At this point, the e-expansion for our model is on no worse formal footing than that for 
the ordinary Ising model. We believe, the asymptotic expansion is valid for both models, 
despite the fact that their critical exponents are different: the exponents for the Ising model 
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in e = 4 — d and the exponents for our model in e = 6 — d are different analytic functions 
with the same asymptotic expansion. The e-expansion can not be used to decide whether 
the lower critical dimension is at e = 3 or at e = 4. 

We conclude that because of instanton corrections the dimensional reduction breaks down 
for the equilibrium RFIM as well as for our nonequilibrium, deterministic zero temperature 
RFIM. In addition there is another "embarrassing" source of non-perturbative corrections in 
the equilibrium RFIM, which we do not have in our problem. There is no reason to expect 
our exponents to be the same as those of the equilibrium RFIM [86] , though the perturbation 
series can be mapped. There might actually be three different underlying functions for the 
same e-expansion for any exponent: one for the pure Ising model, one for the equilibrium 
random field Ising model, and one for our model, so that the exponents in all three models 
would still be different although their e-expansions are the same. 



IX. e-EXPANSION FOR THE AVALANCHE EXPONENTS 

The exponents whose e-expansion we have calculated so far using the mapping to the 
equilibrium RFIM are u, rj, fj, (3, (35, 9 and z. Unfortunately, we cannot extract the avalanche 
exponents r, 1/cr, and 6 from this mapping. The two exponent relations involving these 
exponents 

r - 2 = (7/3(1 - 5) (138) 

and 

l/a = (d-9)v-(3 (139) 

are not enough to determine all three exponents from the information already obtained. 

In the following we will compute r and o directly in an e-expansion. The method 
employed makes use of the scaling of the higher moments of the avalanche size distribution. 
They are being calculated using n replicas of the system for the n'th moment. 
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A. The second moment of the avalanche size distribution 



When calculating 77 and v we have already used all information from the scaling behavior 
of the first moment (S) of the avalanche size distribution: In [129,126] we show that (S) 
scales as the spatial integral over the avalanche-response-correlation function, which in turn 
scales as the "upward susceptibility" dMj dh calculated consistently with the history of the 
system. 

In the Martin-Siggia-Rose formalism it is given by 20 

(S) ~ J dt J d d x(s(t ,x )s(t,x)) f = J dt J d d x(5s(t,x)/5e(t ,x )) f . (141) 

The second moment (S* 2 ) of the avalanche size distribution is the random field average 
of the squared avalanche response. Note that it is not simply the square the expression in 
eq. (141) for the first moment — the product rule for taking derivatives gets in the way: A 
quantity such as 

(s(*o, x )s(t 1 , ari)s(t 2 , x )s(t 3 , x 3 )) f = A + B (142) 



not only contains the term which we need, namely 

&s(ii,zi) Ss(t 3 ,x 3 ) 



5e(t ,x ) 5e(t 2 ,x ] 

but also the terms 



)/ (143) 



20 As we explain in section VIA and in [129,126] the expression 

J dt {s(to,x )s(t,x)) f (140) 

gives the random-field-averaged static response of the system at time t and position x to a positive 
^-function pulse applied at position xq an infinitely long time before t (so that all transients have 
died away). The integral over all space will then give the total static random- field averaged response 
of the system to a ^-function pulse applied at site Xq. This should scale in the same way as the 
first moment of the avalanche size distribution, i.e. as the average avalanche size. 
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S 2 s(t 1 ,x 1 ) 5 2 s(t 3 ,x 3 ) 

B = \TT+ u7^ Mt 3 ,x 3 ))f+ s(ti,a:i)— — — ->/, 144) 

de(t , x )de(t 2 , x ) de(t , x Q )de(t 2 , x ) 



which are not related to (S 2 ). 

How can we separate A from Bl 
One possibility is to introduce a second replica of the system with the identical configuration 
of random fields, the same initial conditions, and the same history of the external magnetic 
field. One can then calculate the response in each of the two replicas separately, multiply 
the results and afterwards take the average over the random fields. Denoting the quantities 
in the first replica with superscript a and those in the second replica with superscript (5 one 
obtains 

A 2 = (s a {t ,x )s a {t 1 ,x 1 )^(t 2 ,x )s p (t 3 ,x 3 )) f 

Ss a (t 1 ,xi) Ss (t 3 ,x 3 ) . ( 
K 5e«{t G ,x )5eP{t 2 ,x ) )f ' 1 } 



since 



and 



¥- = 0. (147) 



Similarly for the n'th moment (S n ) of the avalanche size distribution one would use n 
replicas of the system. In appendix F we make this argument more precise and derive the 
scaling relation between (S 2 ) and A 2 

(S 2 ) f ~ Jdtxj dt a dtpd d x a d d xp(s a (t a ,x )s a (t ,x a )sP(tp,x )sP(t 1 ,x l3 )) f . (148) 

In the following we generalize the RG treatment from previous sections to the case of two 
replicas, and extract the scaling behavior of (S 2 ) from eq. (148) near the critical point. We 
will compare the result to the scaling relation 

(S 2 ) = J S 2 D(S, r, h)dS ~ J S 2 /S T V ± (Sr 1 ^, h/r^dS 

~r( T - 3 V (T S { ± ) (h/r f3S ), (149) 
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where S± is the corresponding scaling function, and obtain the missing information to 
compute the exponents r and a. 



B. Formalism for two replicas 

The generalization of the MSR generating functional to two replicas is rather straight- 
forward. The equation of motion for each spin is the same in both replicas. 

d t s?/r -5H(s a )/5s? = (150) 

and 

dts? /T - 5H( S P)/5 S P = (151) 

where the Hamiltonian TC is given by eq. (32). 

The new generating functional is a double path integral over two 5-functions which 
impose the equations of motion for both replicas. Again we can write the 5-functions in 
their "Fourier" representation by introducing two auxiliary fields s a and s 13 . 

One obtains simply the square of the generating functional from eq. (34), expressed in 
terms of two replicas: 

Z a ? = J J[ds a }[ds a }J J[ds /3 }[ds p ]J[s a }J[s /3 } 

exp(tJ2 J dts«(t)(d t s«(t)/T -5H(s a )/5s«(t))) 
j 

exp (i W dt^(t)(d t s^t)/T - 5H{sP)/5 S P{t))) . (152) 

3 

We note that the two replicas do not interact before the average over the random fields 
is taken. Since Z = 1 we can again average Z directly over the random fields. 

We rewrite the action using the same kinds of transformations to the local fields r] a , fj , 
ff and ff which we introduced previously (see eq. (46)) i.e. 

Z°* = I [dndf} a }[df}[df, P ] 112^,%,$,%] (153) 

J j 
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where Zj[fj", fj'*, fjj, 77^] is a local functional 



= J [ds a }[ds a }[ds^}(expS:l) f , (154) 



and 



~Sf ffj = ~jj dt y^(t)s](t) + is«{t) (d t s«(t) -J^-H-l 

+jjdt{jff j (t) 8 $(t) 



+ 



5V C 



3 5sf 



+ «?(f) \d t s^t) - Jfg - H - fj + ^ I }. (155) 

Here V a and V are given by the linear cusp potential V defined in eq. (31), to be understood 
as a function of s a and s' 3 respectively. 

Again, we expand the action around its stationary point. It is specified by four coupled 
equations, which turn out to be solved self consistently by the replica symmetric mean field 
solution, which we found earlier when studying just one replica: 

C = 0, (156) 



%=0, (157) 



V% = (s a )f, (158) 

fo = (^)f (159) 

Analogously to before [42] we will now expand around the mean-field solution t]q , fj^ 
(JfiQ and Jf/Q denote the local field configurations about which the log of the integrand in 
equation (153) is stationary). Introducing shifted fields r] a = f] a — rfc so that (i] a )f = 0, and 
ff = fj (and correspondingly for if ', and f)P), leaves one with the generating functional 

Z = J[dri a ][dfi a ][d^][df]exp(S afS ) (160) 

with an effective action 
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S aP = -E/ dtJ-'Jrjfit^nt) - £/ dtJ^JrfiWtfit) 
j,i j,i 

+ E X/ I I I I / ^1 ' ' ' dt m + n +p+qU mnpq (ti, ...,t m+n+p+q ) 

j m,n,p,q=0 '"'■"'■P-'J- J 

Vj (tm+n+l) ' ' ' ^/f (^m+ra+p) ? 7j 3 (^m+n+p+l) ' ' ' Vj i^m+n+p+q) ■ (161) 

Here, the are the derivatives of In Zj 13 with respect to the fields fj", 77", rj? and 77^ and 

thus are equal to the local, connected responses and correlations in mean-field theory: 

U *P d d _ d _ 9 

mnpq de"(t m+1 ) de«(t m+n ) deP(t m+n+p+1 ) de? (t m+n+p+q ) 

(s a (t 1 ) ■ ■ ■ s a (t m ) S P(t m+n+l ) ■ ■ ■ s p {t m+n+p )) />J)C . (162) 

As before, local [19] (/) means that we do not vary the local field (r/g )j in the mean-field 
equation 

dts«(t) = J{^) 3 {t) +H + f j - ^ + Je a (t) (163) 

when we perturb the replica a with the infinitesimal force Je a (t) (and correspondingly for 
replica ($). The index c to the average in eq. (162) is a reminder that these are connected 
correlation and response functions. In the same way as we discussed in section VI D the force 
Je a (t) is only allowed to increase with time consistently with the history we have chosen. 
From eq. (162) one sees that u 0npq = if n 7^ 0, u mn0q = if q 7^ 0, and u 0n0q = 0, just as 
we had u 0n = in our earlier calculation for just one replica. 

C. Coarse graining transformation 

The coarse graining transformation is defined in the same way as in the single replica 
case. In the appendix F we give the Feynman rules for loop corrections, and derive the 
canonical dimensions of the various operators in the action. 



62 



D. The scaling of the second moment of the avalanche size distribution 



In appendix F we derive the scaling expression (eq. (148)) for the second moment of the 
avalanche size distribution 

(S 2 ) ~ JdhJ dt a dtpd d x a d d xp(s a (t a ,x )s a (t ,x a )sP(t l3 ,x ) S P(t 1 ,xp)) f . (164) 

In order to find the scaling dimension of (S 2 ) we need to know how 

(s a (t a , x )s a (t , x a )s p {t , x )s (3 (t 1 , x p )) f (165) 

scales under coarse graining. The topology of the diagrams permits no 0(e) loop corrections 
to the corresponding vertex function. Since the anomalous dimensions of the external legs 
(i.e. Greens functions) in the two replicas are also zero at 0(e) it is sufficient to use the 
plain field rescalings to extract the scaling behavior of (s a s /3 s a s 13 ) under coarse graining. As 
shown in appendix F one obtains 

(s a s f3 s a s /3 ) ~ A( 2 (^)- 4 ) (166) 

where A is the cutoff in the momentum shell integrals. 

Inserting this result into eq. (164) along with the canonical dimensions of the various 
times [t] ~ hr z and coordinates [x] ~ A -1 , one obtains 

{S 2 } ~ A- (4+2) (167) 

(Formally including the anomalous dimensions r\ = fj = + 0(e 2 ), one obtains (to first 
order in e) (S 2 ) ~ A~( z +( 2 -v) 2 ) . Similarly, one finds for the higher moments (S n ) ~ 

A-((n-l)z+(2-»7)n)_) 

On the other hand, from eq. (149) we know that (S 2 ) ~ r^ T ~ 3 ^ a S±\h/rP 5 ). If we use 
that r has scaling units A 1 ^, and that r — 2 = a/3(l — 5) (see section IV I, and [129,126]), 
we find by comparison with eq. (167) that 1/er = zv + (2 — rj)u to first order in e. One 
gets the same relation from comparing the dimensions for the n'th moment, which scales as 
(S n ) ~ r ( T -( n+ W a S ( ±\h/rP 5 ). 
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E. Results 



We've seen that 

l/ a = zu+ (2-ri)v = 2 + e/3 + 0(e 2 ). (168) 

If one inserts the result for 1/cr = zv + (2 — r])u = 2 + e/3 + 0(e 2 ) into the relation 
r — 2 = (T/3(l — 5), one obtains 

r = 3/2 + 0(e 2 ). (169) 
From the violated hyperscaling relation 1/a = (d — Q)v — (3 one finds 

Ov = l/2-e/6 + 0(e 2 ). (170) 
This concludes the perturbative approach to the problem. 

X. COMPARISON WITH NUMERICAL SIMULATIONS IN 2, 3, 4, AND 5 

DIMENSIONS 

Figure 15 shows a comparison between the theoretical predictions for various exponents 
and their values as obtained from numerical simulations in 2, 3, 4, and 5 dimensions. These 
prepublication results are courtesy of Olga Perkovic. A complete list of the numerical ex- 
ponents that were measured in the simulations, and a detained description of the algorithm 
that allowed to simulate systems with up to 1000 3 spins is given in a forthcoming publica- 
tion [45]. A quantitative comparison of the results to experiments can be found in [126,128]. 
Some first results and conjectures about the behavior in in two dimensions, which is likely 
the lower critical dimension of our critical point, are presented elsewhere [98]. As is seen 
in the figure, the agreement between the numerics and the results from the e expansion is 
surprisingly good, even down to e = 3. 

The numerical values in 3 dimensions for j3, (35, v and rj seem to have overlapping 
error bars with the corresponding exponents of the equilibrium RFIM. Maritan et al [86] 
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conjectured that the exponents might be equal in a comment to our first publication [41] on 
this system. Why this should be is by no means obvious. The physical states probed by the 
two systems are very different. While the equilibrium RFIM will be in the lowest free energy 
state, our system will be in a history dependent metastable state. Nevertheless, as we have 
seen, the perturbation expansions for the critical exponents can be mapped onto another 
to all orders in e. In appendix G 3 we discuss possible connections between the two models 
that might become clear if temperature fluctuations are introduced in our zero-temperature 
avalanche model. 
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APPENDIX A: HARD-SPIN AND SOFT-SPIN MEAN-FIELD THEORY 



1. Hard-Spin Mean-Field Theory 

In this appendix we derive the scaling forms near the critical point for the magnetization 
and the avalanche size distribution in the hard-spin mean-field theory. At the end we briefly 
discuss changes of nonuniversal quantities for the soft-spin mean-field theory. 

We start from the hard-spin mean-field Hamiltonian: 

H = -J2(JM + H + f t )s z (Al) 

i 

where the interaction with the nearest neighbors from the short range model eq. (3) has been 
replaced by an interaction with the average spin value or magnetization M = (J2 Si)/N. This 
would be the correct Hamiltonian if every spin would interact equally strongly with every 
other spin in the lattice i. e. for infinite range interactions. 

2. Mean-field magnetization curve 

Initially, at H — — oo, all spins are pointing down. The field is slowly increased to some 
finite value H. Each spin Sj flips, when it gains energy by doing so, i.e. when its local 
effective field hi = JM + H + fi changes sign. At any given field H all spins with hi < 
fi < - JM-H will still be pointing down, and all spins with hf > <^> f\ > - JM-H 
will be pointing up. Self-consistency requires that M — J p(f)sidf. Therefore 

[—JM-H roc r-JM-H 

M=(-l) P(f)df+ p(f)df =1-2 p(f)df (A2) 

J-oo J-JM-H J-oo 

is the self-consistency equation from which we can extract the mean-field magnetization as 
a function of H. As in the main text p(f) = exp(-f 2 /2R 2 )/(V2nR) is the distribution of 
random fields. 

For R > \f2~fnJ = R c the solution M(H) of eq. (A2) is analytic at all values of H. 
For R = R c there is a critical magnetic field H C (R C ) = where the magnetization curve 
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M(H) has diverging slope. For R < R c the solution M(H) is unique only for H outside a 
certain interval [H l c (R), H^(R)}. For H in the range between the two "coercive fields" H l c (R) 
and H™(R), the equation has three solutions, two that are stable and one that is unstable. 
Unlike equilibrium systems, which will always occupy the solution with the lowest overall 
free energy, our nonequilibrium (zero temperature) system is forced by the local dynamics 
to stay in the current local energy minimum until it is destabilized by the external magnetic 
field. For increasing external magnetic field this implies that the system will always occupy 
the metastable state with the lowest possible magnetization. Conversely, for decreasing 
external magnetic field, the system will occupy the metastable state with highest possible 
magnetization. One obtains a hysteresis loop for M(H) with a jump, or "infinite avalanche" 
at the upper and lower coercive fields H™(R) and H l c (R) respectively (see figure 5). 
From eq. (A2) follows that dM/dH = 2p(x)/(-t(x)) with x = -JM - H and 

t(x) = 2Jp(-JM - H) - 1 . (A3) 

dM/dH diverges if t(x c ) = 0, which defines x c . For R > R c the slope dM/dH is always 
finite, and t < at all H e [— oo,+oo]. For R < R c however the condition t(x c ) = is 
fulfilled at the critical field H C (R). To obtain potential scaling behavior near this point, we 
expand t(x) around x c = —JM(H C (R)) — H C (R) 

t = 2Jp(x) - 1 = 2J(p(x L ) - 1) + p'(x c )(x - x c ) + l/2p"(x c )(x - x c f + ■■■, (A4) 

where 

2 Jp(x c ) -1 = 0. (A5) 

Then 

X = dM/dH = (-p(x c ))/(J(p'(x c )(x - x c ) + l/2p"(x c )(x - x c f + ■■■)). (A6) 

For a general analytic distribution of random fields p(x) with one maximum with nonvanish- 
ing second derivative (p"(x c ) < 0) this suggests two different scaling behaviors corresponding 
to the cases p'(x c ) = and p'(x c ) ^ 0. 

67 



Let us consider the case p'(x c ) = first. For a Gaussian distribution of width R = R C with 
zero mean this implies that x c = —JM(H C ) — H c = and consequently p(x c ) = 1/ (\/2ttR c ) 



and p"(x c ) = l/(2irRl). With eq. (A5) one obtains R c = This is in fact the largest 

possible value of R for which M(H) has a point of diverging slope. 

Integrating eq. (A6) leads to a cubic equation for M and the leading order scaling be- 
havior 

M{r, h) ~ \rfM±{h/\rf 5 ) , (A7) 

for small h = H — H(R C ) and r = (R c — R)/R. In mean-field theory 5 = 3 and A4± is given 
by the smallest real root g±(y) of the cubic equation 

o 12 12y/2 /A . 

' T 7'-^» = 9 ' (A8) 

where here and throughout ± refers to the sign of r. 

The other case (p(x c ) = 1/(2 J) and p'(x c ) ^ 0) is found for distributions with R < R c . 
Integrating eq. (A6) with x c = —JM(H C (R)) — H C (R) yields a quadratic equation for the 
magnetization and the scaling behavior 

M - M(H C (R)) ~(H- H C (R))< (A9) 

with C = 1/2 for H close to H C {R). From eq. (A2) and eq. (A5) one finds H C (R C ) = 0, 
H C (R) ~ r 135 for small r > 0, with (35 = 3/2 and H C (R = 0) = J. The corresponding phase 
diagram was shown in figure 7. 

Note that the scaling results for R close to R c as given in eq. (A7) remind one of the 
scaling results of the Curie- Weiss mean-field theory for the equilibrium Ising model near the 
Curie-temperature (T = T c ). For T < T c however, the equilibrium model has a discontinuity 
in the magnetization at H = 0, while for R < R c our model displays a jump in the mag- 
netization at a (history dependent) nonzero magnetic field H C (R), where the corresponding 
metastable solution becomes unstable. Our infinite avalanche line H C (R) is in fact similar 
to the spinodal line in spinodal decomposition [53]. 
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Note also that this mean-field theory does not show any hysteresis for R > R c (see figure 
5). This is only an artifact of its particularly simple structure and not a universal feature. 
For example, the continuous spin model, which we use for the RG description in section VI 
has the same exponents in mean-field theory, but shows hysteresis for all disorders R, even 
for R > R c . In that case there are two critical fields H^(R C ) and H l c (R c ) — one for each 
branch of the hysteresis curve (see figure 6, and figure 8). 

3. Mean-field avalanche-size distribution 

As we have already discussed in the main text, one finds avalanches of spin flips as the 
external field is raised adiabatically. Due to the ferromagnetic interaction a flipping spin 
may cause some of its nearest neighbors to flip also, which may in turn trigger some of their 
neighbors, and so on. In mean-field theory, where all spins act as nearest neighbors with 
coupling J/N, a spin flip from —1 to +1 changes the effective field of all other spins by 2J/N. 
For large N, the average number of secondary spins that will be triggered to flip in response 
to this change in the effective local field is then given by n tr ig = (2J/N)Np(—JM — H) — 
2Jp(—JM — H). If n t rig < 1, any avalanche will eventually peter out, and even in an infinite 
system all avalanches will only be of finite size. If n trig = 1, the avalanche will be able to 
sweep the whole system, since each flipping spin triggers on average one other spin. This 
happens when the magnetic field H takes a value at the infinite avalanche line H = H C (R), 
with R < R c , since n tr ig = 1 is equivalent to t — (see eq. (A3)). 

Considering all possible configurations of random fields, there is a probability distribution 
for the number S of spins that flip in an avalanche. It can be estimated for avalanches in 
large systems, i.e. for S « N: For an avalanche of size S to happen, given that the primary 
spin has random field f, it is necessary that there are exactly S — 1 secondary spins with 
their random fields in the interval [fi, fi + 2(J/N)S]. Assuming that the probability density 
of random fields is approximately constant over this interval, the probability P(S) for a 
corresponding configuration of random fields is given by the Poisson distribution, with the 
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average value A = 2JSp(-JM - H) = S(t + 1), where t = 2Jp(-JM - H) - 1: 

P(S) = exp(-A). (A10) 

This includes cases in which the random fields of the s spins are arranged in the interval 
[ft, ft + 25'( J/N)] in such a way that they do not flip in one big avalanche, but rather in 
two separate avalanches triggered at slightly different external magnetic fields. Imposing 
periodic boundary conditions on the interval [f i: fi + 2S(J/N)] one can see that for any 
arrangement of the random fields in the interval there is exactly one spin which can trigger 
the rest in one big avalanche. In 1/S of the cases, the random field of this particular spin 
to trigger the avalanche will be the one with the lowest random field, as desired. Therefore 
we need to multiply P(S) by 1/S to obtain the probability D(S,t) for an avalanche of size 
S starting with a spin flip at random field fi = —JM — H 

D(S, t) = S is - 2) /(S - l)!(t + l)( 5 - 1 ) e - s (*+ 1 ) . (All) 

With Stirling's formula we find for large S the scaling form 

D(S, t) ~ -=L^ exp(-^ 2 /2) • (A12) 



To obtain the scaling behavior near the two different critical points, we will insert into the 
expression in eq. (A12) the expansion of t(x) around x c from eq. (A4). 



4. Avalanches near the critical endpoint 

Near the critical point (R c , H C (R C )), where x c = and p'(x c ) = we obtain t = 2 J(p(0) — 
1) + Jp"(0)(—JM — H) 2 , which (by equation eq. (A7)) leads to the scaling form 

t^ r [lTl/^g±(h/\rf s ) 2 ]. (A13) 

g± was defined in eq. (A8), and ± again refers to the sign of r = (R c — R)/R. Inserting the 
result into eq. (A12) one obtains 
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D{S,r,h)~S- r V±(Sl\r\- 1 l°,hl\r\ f,S ), (A14) 
with the mean-field results t = 3/2, a = 1/2, (35 = 3/2, and the mean- field scaling function 

V ± {x, y) = -M^fe) 2 ] 2 /? . (A15) 

V 27T 



5. Mean-field avalanche size distribution near the oo-avalanche line ("spinodal line") 

For R < R c one has p'(x c ) ^ 0, so that the expansion for t becomes 

t = 2J{p{x c ) - 1) + 2Jp'{x c ){x - Xc ) + ... = 2Jp'{x c ){x - x c ) (A16) 
= 2Jp'(x c )(-J(M - M{H C {R))) - (if - H C (R))) 

Following the steps that led to eq. (A9) we arrive at 

t = -2y/jf/(x c )(H - H C (R)) + higher orders m (H - H C (R)) (A17) 

so that for H close to the onset to infinite avalanche (with H < H™(R) for increasing field 
H and H > H l c (R) for decreasing field) 

D(S, (H - H C (R))) ~ -=L^ exp {-2[p'{- JM - H)J\S\H - H C {R)\) . (A18) 



or, written more generally, 

D{S,H - H C (R)) ~ \/S T F{S\H - H C {R)\^ K ) (A19) 

with k — 1 and r = 3/2 in mean-field theory, and T the corresponding mean-field scaling 
function. 

6. Modifications for the soft-spin mean-field theory 

In section VI we have for calculational convenience switched from the hard-spin model, 
where each spin Sj could only take the values ±1, to a soft-spin model, where s$ can take 
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any value between — oo and +00. In realistic systems these soft-spins can be considered as 
coarse grained versions of the elementary spins. The corresponding Hamiltonian with the 
newly introduced double-well potential 



V(si) = 



(A20) 



k/2 ( Si + l) 2 for s < 
fc/2 (s { - l) 2 for s > 

to mimic the two spin states of the hard-spin model, was given in eq. (32). In the mean- 
field approximation, where the coupling term —JijSiSj is replaced by — J2i JMsi with M = 
J2j s j/N, we obtain 



H = - £ {( JM + H + ftsi - V( Si )} . 



(A21) 



For adiabatically increasing external magnetic field the local dynamics introduced earlier 
implies that each spin will be negative so long as the lower well Hamiltonian 



H- = fc/2 ( Si + l) 2 -(H + JM + fi) Si 



(A22) 



does have a local minimum with 5H/5s = for negative Sj. This implies that < if 
5 



Ssi 



[k/2 { Si + I) 2 -(H + fi + JM) Si } s=0 > , 



(A23) 



else Si will be stable only at the bottom of the positive potential well, where 
-L-H+ = /-[fc/2 ( Si - l) 2 -(H + JM + fi) Si } = . 

OSi OSi 

We conclude that for the given history 

Si < for ft < - JM -H + k 
Si>0 for fi > - JM -H + k. 

From the self-consistency condition 



(A24) 



(A25) 



( Si ) EE J P (fi) Si dfi = M 



(A26) 



we derive the mean-field self-consistency equation for the soft-spin magnetization M u (for 
increasing external magnetic field): 
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r-JM-H+k 

M U (H) = (k + H)/(k -J)- 2k/ (k - J) / p(f)df . (A27) 

J — oo 

Correspondingly one finds for the branch of decreasing external magnetic field: 

r -JM-H-k 

M^H) = (k + H)/(k - J) - 2k/(k - J) / p(f)df . (A28) 

J— oo 

Figure 6 shows the corresponding hysteresis loops in the three disorder regimes R < R c = 
ix.J{k/{k — J)), where the hysteresis loop has a jump, R = R c , where the jump has 
shrunk to a single point of infinite slope dM/dH, and R > R C) where the hysteresis loop is 
smooth. In contrast to the hard spin model, this model displays hysteresis even for R > R c . 

The critical magnetic fields H^(R) and H l c (R) at which the slope of the static magne- 
tization curve diverges are found by differentiating eqs. (A27) and (A28) with respect to 
H and by solving for dM sta t/dH. One finds (for increasing external magnetic field) that 
dM s tat/dH ~ 1/r with 

r = (2k/ (k - J))Jp(-JM stat -H + k)-l. (A29) 

r is defined analogously to the paramter t in eq. (A3). (It is worth mentioning that t = 
X~ 1 k/(J(k — J)), where x _1 is the constant term in the propagator of the RG treatment, 
which was introduced in eq. (73) in the main text.) The critical field H C (R) is given by 
the solution to the condition r = 0. To find the scaling behavior near the critical point 
one can expand eq. (A27) around H^(R), and correspondingly eq. (A28) around H l c (R). 
For increasing external magnetic field the critical point R = R C) H = H^(R C ) and M = 
M c = M U (H^(R C )) is characterized by the equation r = and p'(—JM c — H c + k) = 0, i.e. 
-JM C - H c + k — 0. It follows that R c = -^fzj- Inserting these results into eq. (A27) one 
obtains M" = 1 and H^(R C ) = k — J. Similarly for a decreasing external magnetic field one 
finds H l c (R c ) = —(k — J) and M l c = Mi(H l c (R c )) = —1 The corresponding modified phase 
diagram is depicted in figure 8, with H™(R = 0) = +k and H l c (R = 0) = —k. 

In the same way as discussed in appendix A for the (static) hard-spin model, expanding 
eqs. (A27) and (A28) around M c , H c and R c yields cubic equations for the magnetization 
and one obtains the scaling behavior near the critical point 
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m ~ /i 1/3 (A30) 

for m = M — M c and h = H — H C (R C ). Consequently the slope of the magnetization as a 
function of h diverges at the critical point 

dm/dh ~ /T 2/3 . (A31) 

It turns out that in fact none of the universal scaling features we discussed for the hard 
spin model is changed. The mean-field critical exponents f3, 5, r and a, and the scaling forms 
near the critical point are the same as in the hard-spin model. (A "spin-flip" in the hard 
spin model corresponds to a spin moving from the lower to the upper potential well in the 
soft-spin model.) Apart from modifying some nonuniversal constants, the new parameter 
k > J in the definition of the soft-spin potential does not appear to change the calculation 
in any important way. 

a. Soft-Spin Mean-Field Theory at Finite Sweeping Frequency Q 

In section VIA, eq. (62), we have derived the following equation of motion for each spin 
in the dynamical soft-spin mean-field theory, as the external magnetic field H(t) = H + flt 
is slowly increased 

±-d tSj (t) = J V %t) +H + fj- ^ + Je(t) . (A32) 

With the definition of the potential V from eq. (31) this becomes 

l/(T k)d tSj (t) = - Sj (t) + JVj(t)/k + H(t)/k + fj/k + sgn(s 3 ) + Je{t)/k . (A33) 

From eq. (54) we know that i]°(t) = (s) = M(t) is the time dependent mean-field magneti- 
zation of the system. It can be calculated by taking the random-field average of eq. (A33) 
and solving the resulting equation of motion for rj°(t). One can show [126] that for driving 
rate fl/k small compared to the relaxation rate kTo of the system, for all values of Hq the 
solution i] (t) can be expanded in terms of (fl/T ) in the form 
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V °(t) = M(t) = M stat (H ) + (fi/r o r/i(#o)t + (n/To) P2 f 2 (Ho)t 2 + ■■■ (A34) 

with < pi < P2 < ■ ■ •■ The Pi depend on whether R < R c or R = R c . M stat (H Q ) is 
the solution of the static mean-field theory equation (A27) for the given history. If the 
series converges for Q — > 0, it follows that r]°(t) approaches the constant magnetization 
M sta t(H ) in the adiabatic limit. This is certainly expected for H Q away from the critical 
field H C (R), where the static magnetization is non-singular: as Q tends to zero the time 
dependent magnetization M(t) simply lags less and less behind the static value M stat (H(t)). 
The magnetization M(t) can be expanded as M(t) = M stat (H ) + [dM/dH] Ho Qt + ■ ■ ■ and 
converges towards M stat (H ) as f2 — > 0, as long as all derivatives [d n M stat /dH n ] Ho are well 
defined and finite. This argument however does not obviously hold at the critical fields 
H = H C (R) with R < R c , where dM sta t/dH and all higher derivatives diverge. Using 
boundary layer theory one can show [126] that even at these singular points M{t) converges 
toward its static limit M(H C (R)) as f2 — > 0, though with power laws smaller than one in 
Q, as indicated in eq. (A34). Since we use M stat (H ) as the foundation for our e-expansion, 
this is reassuring. 

APPENDIX B: TILTING OF THE SCALING AXES 

In the appendix on mean-field theory we derived scaling forms for the magnetization 
and the avalanche size distribution, which depended on the (mean-field) scaling fields r = 
(R — R c ) I R and h = H — H c . In finite dimensions however, the corresponding scaling forms 
may depend not on r and h, but on rotated variables 

r' = r + ah (Bl) 

and 

ti = h + br. (B2) 

The amount by which the scaling axes r' — and h! = are tilted relative to the (r, 0) and 
(0, h) direction in the (r, h) plane is a nonuniversal quantity and has no effect on the critical 
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exponents. Nevertheless it can be important in the data analysis. (The numerical results in 
3, 4, and 5 dimensions do indeed seem to indicate a slight tilting [45].) 

1. Extracting (3 and 06 from the magnetization curves 

With m(r', h!) = M(r', h!) - M c , where M c = M(0, 0) is the magnetization at the critical 
point, we obtain the scaling form 

m~ (r'fM±(ti/(r'f s ). (B3) 

In simulations and experiments however the magnetization is given as a function of r and h 
rather than r' and h'. We rewrite eq. (B3) in terms of r and h by inserting the definitions 
of r' and h! from eq. (Bl) and eq. (B2) and use the fact that (35 > 1 in the cases we are 
considering. One obtains for the leading order scaling behavior 

M - (M c + cr) ~ r l3 M±({h + br)/r^ s ) (B4) 

(where c is a constant and Ai± is the appropriate function). Equivalently this can be written 

as 

dM/dH ~r p - ps M±{{h + br)/r f)S ). (B5) 

To extract the critical exponents f3 and (56 one can then use collapses of dMj dH in the 
same way as if the mean-field scaling forms were valid in finite dimensions, except for the 
presence of a new tilting parameter b, that has to be varied to its correct value simultaneously 
with (3 and (38 to find the best collapse of the data curves. (A more detailed description of 
the procedure is given in reference [45]). 

2. Extracting the correlation length exponents v and u/{(35) 

Similarly to eq. (B3) the scaling for the correlation length in finite dimensions takes the 
form 
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t,~{r')- v y±{til{rT), (B6) 

Thus £ ~ (h')~ v ^^ along r' — and £ ~ (r') - " along ft/ = 0. Figure 16 shows contour 
lines in the (r',h') plane for the correlation length. Since v/{f3$) < v, £ changes faster in 
the (0, h!) direction than in the (r', 0) direction. This implies that the correlation length 
diverges with the dominant exponent v/((38) when the critical point is approached along any 
direction other than hi = 0. This can be used to extract v/{f38) from collapses of numerical 
or experimental curves for the correlation funcion G(x, h, r) measured as a function of h at 
fixed R = R c : 

G(x, h, r = 0) ~ l/x^- 2+r '^g ± (xh-^^ 5 \0) (B7) 

with the appropriate scaling function Q±. (Even in 3 dimensions v/(/36) < v is still correct 
and the tilting is small.) 

On the other hand it seems rather difficult to find the weaker scaling direction (r', 0) 
accurately enough to be able to extract v by approaching the critical point along this line. 
If, instead, we integrate the correlation function or the avalanche size distribution over 
h', we obtain a scaling form that depends only on r'. From collapses of the integrated 
functions, the exponent v can then be extracted, even without knowledge of the exact size 
of the tilting angle. Practically, in the analysis of our simulation data we actually integrate 
over the magnetic field H rather than h'. On long length scales this is equivalent to an 
integration over h': The integration path r = const = c is written in the (h',r f ) plane 
through r' = a±c + b\h and h' = a^c + b 2 h. Under coarse graining the h! component grows 
faster than the r' component, so that after several coarse graining steps the integration 
path is deformed into a straight line parallel to h! . This procedure yields the dominant 
contribution to an integrated scaling function G int (x,r) = J dhG(x,r, h) in the region near 
the critical point. It allows us to treat r as a scaling field and to extract the weaker exponent 
v. we have £ ~ r~ v after integration over the hysteresis loop. For the integrated avalanche 
correlation function, for example, this means G int (x,r) ~ l/x d+ M v gg t (xr- u ), with the 
appropriate scaling function Q™*. 
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The rotation of the scaling axes (see also [45]) is not apparent from the e-expansion. For 
the renormalization group calculation we have linearized the coarse graining transformation 
around the fixed point. The rotation of the scaling axes is due to to nonlinear corrections 
introduced during the flow to the fixed point. 

From the mapping of the perturbation series for the critical exponents in our model 
to the perturbation series for the critical exponents in the equilibrium random field Ising 
model, one might expect that for R > R c the H — > (-H) symmetry of the fixed point in the 
thermal random field Ising model is mapped to an {h!) — > {—hi) symmetry in our problem. 
(Nonperturbative corrections might destroy this of course.) Note that this symmetry only 
emerges on long length scales, near the critical point, where the mapping to the thermal 
random field Ising model holds. 

APPENDIX C: SOME DETAILS OF THE RG CALCULATION 

1. Calculating some u mn coefficients 

In section VIA in eq. (61) we have given an expression for the coefficients u mn in the 
expansion around mean-field theory: 

d d 



u. 



((s(fx) - V °{ti)) ■ ■ ■ {s{t m ) - ^(U))^ , (CI) 



de{t m+ i) de{t m+n 
where Sj{t) is the solution of the local mean-field equation 

l/{r k)d tSj {t) = - Sj (t) + Jrf^/k + H{t)/k + fj/k + sgn{ S j) + Je{t)/k . (C2) 

We need to insert the solution for rj°{t) from appendix A 6 into eq. (C2) to calculate the 
higher response and correlation functions u mn as given in eq. (CI). 



a. The coefficient uifi 

The vertex u lj0 is defined as 
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ui fl (t) = Mt))i-V°(t). (C3) 

From the stationary phase equation eq. (54) we have rj°(t) = (sj(i))j. Therefore (by con- 
struction) 

«!,(>(*) = 0- (C4) 

b. The coefficient 1/1,1 (ii, t 2 ) 

As is shown in appendix A 6, rf{t) can be expanded in terms of Q (at least for R > R c 
and for R < R c before the jump up to H C (R))\ 

rf{t) = M(H ) + il p t + ---, (C5) 

where M(H ) is the static magnetization, p > 0, and • • • implies higher orders in Q. In- 
serting this expansion into eq. (C2) and eq. (CI) allows us to calculate the coefficients u mn 
perturbatively in fl. Only the lowest order remains as Q — > 0. The vertex function Un(t 1 , t 2 ) 
is then given by 

^mn^o[9t(Hm e _ +0 (s(i2)|if(t 2 )+jee(t 2 -ti) - s{t 2 )\ H (t 2 )) f /^)] ■ (C6) 

To evaluate s(t2)\H(t 2 )+Jee(t 2 ~t 1 ) from the equation of motion eq. (C2) we have to consider 
three cases: 

1. Neither Sj(t)\H(t) nor s(t)\H(t)+Je9(t-ti) a ^ an y time t with ti < t < t 2 . To lowest 
order in Q the solution of eq. (C2) is then given by simple relaxation: 

sY flip {t)\ H{t)+ j ee{ t- tl) - s 3 (t)\H(t) = Wk)(l - exp [-kV (t - h)]^ - t ± ) . (C7) 

2. The unperturbed Sj(t)\jf(t) AiP s a ^ some time tj with t\ < tj < t 2 . The fraction 
of spins for which this is the case is proportional to Q. Thus, for Q — > these spins will 
yield no contribution to Uu for finite t 2 — ti.( This is true even if one chooses to keep 
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H(t 2 ) — H{t\) = AH 7^ fixed as f2 — > 0. One finds that the resulting contribution to tin 
involves terms of the form exp[— kT AH/Q] which clearly vanish as f2 — > 0.) 

3. The unperturbed spin (*) | J? (t) does not flip at any time between t\ and t 2 , but the 
perturbed spin s(t)\H(t)+Je6(t-ti) does flip at time tj with t\ < tj < t 2 . For fixed, finite 
t 2 — ti we can again expand the contribution As to ti u in terms of Q. In the adiabatic limit 
at fixed, finite t 2 — t\ only the lowest order term survives, which is independent of Q. It 
is calculated at constant rj°(t 2 ) = rj°(ti) and H(t 2 ) = H(ti) = H, as done in the following 
paragraph: 

The time t j with t\ < t j < t 2 at which a spin moves from the lower to the upper well is 
given by 

s(tj)\ H+ jee(t J -t 1 ) = 0. (C8) 

Before the perturbation is switched on, the static solution for the spin in the lower potential 
well is 

s i (t 1 ) = Jr ] /k + H/k + f i /k-l. (C9) 

As soon as the perturbation is switched on, the spin starts relaxing into its new equilibrium 
position according to eq. (C7). Its value will remain negative until time tj, when it will flip 
to a positive value in the upper potential well. Thus (using eq. (C9) and eq. (C7)) 

= Si (tj) = (1 - exp(-fc(*j - h)))Je/k + {Jrf/k + H/k + fc/k - 1) , (CIO) 

or 

t J = t l - l/ATo ln(l + (Jrf/k + H/k + fc/k - l)k/Je) . (Cll) 

The shift in sj(t) at time t 2 > tj will therefore not only consist of a contribution due 
to simple relaxation as given by eq. (C7) but also a contribution due to the spin flip. It is 
proportional to the distance of the equilibrium point of the upper well to the equilibrium 
point of the lower well, which is at all times two in our model. Solving eq. (C2) with 
sgn(si(t < tj)) = —1 and sgn(si(t > tj)) — +1, we find in this case 
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s{ Hp (t 2 )\ H+Mt2) - Si (t 2 )\ H = Je/k(l - exp(-kF (t 2 - ii)))6(t 2 - *i) 
+2(1 - exp(-kF (t 2 - tj)))e(t 2 - tj)e(tj - h) . (C12) 

The first of the two terms on the right hand side of the equation also appears in eq. (C7). 
We can therefore separate the contribution solely due to the spin flip from the relaxational 
part: 

From eq. (C7) and eq. (C12) one obtains 

As = s{ lip (t 2 ) - sT fUp {t 2 ) = 2(1 - exp(-AT (i 2 - tj))) . (C13) 
The disorder averaged shift in the local magnetization is then 

(s(t 2 )\ H+Mt) - s(t 2 )\ H ) f ~ (sr fHp (h)\ H+<t) - Sl (t 2 )\ H ) (C14) 

+ f 2 dtjdN{tj)/dtj2{l - exp{-kT {t 2 - tj)))Q{t 2 - tj) . 
Jti 

dN(tj) is the number of spins that are going to flip in the time interval [tj,tj + dtj] with 
ti < tj < t 2 . From eq. (CIO) we find a relation between the random field and tf. 

/. = -Jrf -H + k- Je(l - exp(-kT (tj - h))) . (C15) 

Therefore 

dN{tj) = (C16) 
dtj\dfi/dtMfi = -Jr]° -H + k- Je(l - exp(-A;r (t J - h)))) 
= dtj(Jekr exp(—kr (tj — ti)) 

p(fi = -Jrf -H + k- Je(l - exp(-ATo(tj - h)))) 
= dtj(JekT exp -kF (tj - tMfr = -Jrf -H + k)) + 0{t 2 ). 

Eq. (C7), eq. (C12), and eq. (C14) then yield for the average shift (for t 2 > ti) divided 
by e, in the limit of small e 

lim e _>o[(si(t 2 )\H+je(t2) ~ Si{t 2 )\ H )\/e (C17) 
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= [J/k + 2Jp(f l = -Jrf -H + k) + 

exp(-AT (t 2 - h)){-J/k - 2J(1 + kT (t 2 - fi)) 
p(fi = -Jv°-H + k))]G(t 2 -t 1 ), 

which correctly goes to zero for t 2 — > ti. (The terms proportional to p are due to spin flips 
and the others stem from simple relaxation.) For {t 2 — t x ) — > oo it approaches the static 
value J/k + 2Jp(fi = Jrf — H + k), which is the static response to a small perturbation. 
The ui t i(ti,t 2 ) vertex is given by the negative time derivative of eq. (C17) with respect to 
h, i.e. 

«i,i(*i> h) = (J/k + 2Jp(f t = -J V ° -H + k))5(t 2 - (C18) 
+5{t 2 - tJlexpi-kFoih - t^i-J/k - 2J(1 + A;r (t 2 - h)) 

p(fi = -Jv°-H + k))] 
+&(t 2 - h)d tl [exp(-kT (t 2 - h))(-J/k - 2J(1 + A;r (t 2 - h)) 

p(fi = -Jv°-H + k))]. 

In the corresponding term in the action this is multiplied by r](t 2 ) and r)(ti) and integrated 
over dti and dt 2 . The terms multiplied by 5(t 2 — ti) in eq. (C18) cancel exactly, so that there 
is only left the term multiplied by 9(t 2 — ti). We perform an integration by parts in t\ and 
obtain two terms: The boundary term which has a purely static integrand and leaves only 
one time integral, and a time dependent part with two time integrals. 
The static term contributing to the action is then 

/+oo 
dhrfrMhX-J/h - 2Jp(f t = -Jrj° - H + k)) . (C19) 
-oo 

The dynamical part can be written as 

/+oo r+oo 
dt 2 / dh&ih - tJrjfoXdttfih)) exp(-kT (t 2 - fx)) 

{-J/k - 2J(1 + AT (t 2 - h)))p{-Jrf -H + k). (C20) 

In the low frequency approximation this becomes 
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/CO 
dhftitJdtMh) l-J/k - Upi-Jr] -H + k)]/(T k) 
-oo 

/oo 
dttfitJdtrttMTo , (C21) 
-oo 

with 

a=[J/k + 4Jp(-Jr]° - H + k)]/k, (C22) 
as can be seen using that 

e -*ro(t2-ti) e ( i2 _ fl ) = A_ J due^-^/ikTo + iuj) , (C23) 

and 

/IcV 
dwe^-^ -—^-. (C24) 

Eq. (C21) contributes to the u iu" term in the propagator expressed in frequency space. 
We have performed the above calculation for the case H(ti) = H(t 2 ). If instead one 
keeps H(t 2 ) — H(t±) = AH fixed as f2 — > (i.e. t 2 — h — > oo), one obtains 

Zim e _ (si(t2)|H+ e (t 2 ) -Si(* 2 )|H)//e= J/k + 2Jp(-Jrj°(t 2 ) - H + k) (C25) 

up to dynamical corrections of the form (exp(— Ai?r /Jl)), which are negligible as f2 — > 0. 
Consequently, the derivative with respect to (— £i) yields zero in this limit. Therefore there 
is no contribution to the action from these cases and the result converges to the expressions 
in eq. (C19) and eq. (C20) as Q -> 0. 



c. The coefficients u± t2 , ui,3 and u 2 $ 

The coefficients u ij2 and u lj3 at field if are calculated similarly. One obtains for the 
terms in the action corresponding to u ij2 in the adiabatic limit: 

/ d d x I dtfj(x,t) (w[i](x,t)} 2 + ( dt 2 a(t,t 2 ,t 2 )r)(x,t 2 )d t2 r)(x,t 2 ) 

J J — CO \ J — CO 

+ I dt 2 [ dtia(t,t 1 ,t 2 )d ta ri(x,t2)dt l ri(x,t 1 )) . (C26) 

J— CO J —CO J 
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Here, w = —2J 2 p'(fi = —Jr/ — H + k), and a(t,ti,t 2 ) is a transient function due to the 
relaxational dynamics of the system. It consists of terms proportional to exp{— T (t — ti)} 
or exp{— T (t — t 2 )}. The transient terms proportional to a(t, ti, t 2 ) turn out to be irrelevant 
for the critical behavior observed on long length scales. 

The static and the transient terms in the action contributed by -Ui )3 , are calculated 
similarly. Again, only the static part turns out to be relevant for the calculation of the 
exponents below the upper critical dimension. It is given by 



f d d x f °° dtufi{x,t)[t]{x,t)} 3 , (C27) 

with u = 2J 3 p"(fi = -J77 -H + k). 

Finally, the vertex U2fi(ti,t 2 ) = (si(ti)si(t2)) is calculated similarly. It is a local corre- 
lation function instead of a local response function. Therefore the times t\ and i 2 can be 
infinitely far apart, i.e. even for H{t\) 7^ Hfa) the vertex u 2t o is still nonzero. One obtains: 



., . / f-H(t 2 )-rf\t 2 )+k 



u 2 fi{t u t 2 ) = R /k 2 +yj p(h)dh 

-H(ti)-ri (ti)+k 



Hz 



p(h)dh 

H(t 2 )-ri (t2)+k , 



( r-H(t 2 )-rP(t2)+k \ 2 ( i—H(t 2 )-rf(t 2 )+k 

-4iy p(h)dh) -4iy (h/k)p{h)dh) 



-2 / (h k)p(h)dh) , (C28) 

\J-H(t 2 )-rP(t 2 )+k ) 

which is positive (or zero) for any normalized distribution p(f). 



2. Feynman Rules 

In the following discussion we denote with u mn the static part of u mn (ti...t m+n ), i.e. the 
part which, (for n ^ 0, after taking the time derivative and integrating by parts) is not 
multiplied by any time derivative of the fields. This is usually the only part of the vertex 
which is not irrelevant under coarse graining (except for the propagator term, which also 
has a contribution proportional to iuj). We formulate the Feynman rules in general terms, 
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although for our discussion only = x 1 / J + 1) u i,2 = w, = u and 1^2,0 come into 
play. 

We introduce the following notation for the Feynman diagrams, which we will use in 
the perturbation expansion : A vertex u^ n is a dot with m outgoing arrows (one for each 
fj operator) and n incoming arrows (one for each r\ operator.) Figure 17 (a) shows the 
graph for u as an example. The propagator can only connect fj(t) and rj(t') operators. It is 
important to remember that causality must be obeyed, i.e. t' > t. Figure 17 (b) shows an 
example of a diagram which violates causality. 

The vertex -u 2 ,o is denoted as shown in figure 17 (c). The black ellipse connects the two 
parts of the vertex that are taken at different times. 

The symmetry factors for each diagram are obtained in the usual way [105,77] by drawing 
all topologically distinct graphs and counting their multiplicities. One must not forget to 
include the factors l/(m\n\) that appear in the expansion (see eq. (60)). 

In each diagram, momentum conservation requires that vertices should be connected by 
loops rather than a single, dead end propagator line. Figure 17 (d) shows an example of a 
diagram that is zero [78]. 

Each internal loop contributes an integral over the internal momentum-shell [53] (see 
also below) 



where A is the cutoff and b > 1. 

Integration over time is already performed in eq. (C29). (In the low frequency approx- 
imation the propagator can be approximately taken to be S(t — t'), when we integrate out 
modes in the infinitesimal momentum shell A/b < q < A [19].) 



As we have mentioned in section VI C it turns out that on long length scales different 
magnetic fields decouple and the static critical exponents can be extracted from a renor- 




(C29) 



3. Implementation of the history 
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malization group analysis performed at a single, fixed value of the external magnetic field 
H Q due to a separation of time scales. In the following paragraph we will show that this 
statement is self-consistent using an argument by Narayan and Middleton in the context of 
the CD W dep inning transition [20] . 

An expansion around mean-field theory in the way performed here corresponds to first 
increasing the magnetic field H within an infinite ranged model and then tuning the elastic 
coupling to a short ranged form, while the actual physical behavior corresponds to first 
tuning the elastic coupling to a short ranged form and then increasing the force within the 
short ranged model. The concern is that in the presence of many metastable states the 
critical behavior of the two approaches might not be the same. For example, spins might 
tend to flip backwards upon reduction of the interaction range in the expansion around 
mean-field theory. Although there will of course always be some spins for which this is the 
case, no such effects are expected on long length scales since the susceptibility is actually 
more divergent near the critical point for d < 6 than in mean-field theory: 

(dm/dh) h=0 ~ r~ 7 (C30) 

and 

{dm/dh) r=0 ~ h 1 ' 5 ' 1 (C31) 

with 

j (in d < 6) > j (in mean field theory) (C32) 

and 

(1 - 1/5) (in d < 6) > (1 - 1/5) (in mean field theory) (C33) 

as we learn from numerics and analytics. This suggests that on long length scales spins 
tend to flip forward rather than backward upon reduction of the coupling range. One would 
therefore expect the expansion around mean-field theory performed here to correctly describe 
the critical behavior. 
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Reassured by this self-consistency argument we now briefly discuss for separated time 
scales that the decoupling of the different magnetic fields is consistent also within the RG 
description. As discussed in appendix C 1 the response functions Mi, n (ii, • • •, t n +i) with fixed 
H(ti) 7^ H(tj) with j'^1 tend to zero in the adiabatic limit. However, the original action S 
of eq. (60) also contains terms of the form u 2>0 (Hi, H 2 ) which do couple different fields even 
as Q — > 0. These "multi- field" vertices however do not contribute to the renormalization of 
the vertices evaluated at a single value of the external magnetic field, because the propagator 
does not couple different field values. It turns out that the multi-field vertices are also irrel- 
evant on long length scales: setting up the RG in section VII B we have treated u 2t o(H , H ) 
as a marginal variable. In fact, we have been choosing the rescaling 7^ of the fields fj in 
rj(x, t) = b~ d l 2 ~ z ~"'^ 2 f] '(x, t) such that u 2 ,o(H , H ) remains marginal to all orders in pertur- 
bation theory. Similarly, the rescaling 7^ of the fields 77 in r)(x,t) = b~ d l 2+2 ~ 1 ' t 'l 2 r]'(x,t) is 
chosen such that the coefficient of the q 2 term remains marginal to all orders also. (These 
choices are made so that the rescaling of the response and the cluster correlation function un- 
der coarse graining immediately gives their respective power law dependence on momentum 
(see section VII A)). To decide what this implies for u 2j0 (Hi, H 2 ) at two different magnetic 
fields Hi 7^ H 2 , we need consider two loop order corrections for u 2 p(Hi, H 2 ). In figure 18 
we have indicated the magnetic fields corresponding to the times at which the vertices are 
evaluated, taking into account that the propagators do not couple different fields. One finds 
the following recursion relation: 

\u' 2fi {H ll H 2 ) = \u 2fi {H ll H 2 ) (C34) 
+ ^u 2fi (H 1 ,H 2 )(u 2fi (H 1 ,H 2 ) 2 -u 2 , (H 1 )u 2fi (H 2 )) 
^ir 6 2 (47r) 6 e 2 /( M2 ,o(i/i)^o(// 2 ))ln6, 

with u 2j0 (Hi) = u 2> o(H 1 , Hi). The term that is subtracted stems from the second order 
correction to the rescaling of the fields r) (from terms as depicted in figure 18 but with 
Hi=H 2 ). 

From eq. (C28) it follows that u 2fi (Hi, H 2 ) > 0. Furthermore, using the same equa- 
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tion one finds after some algebra that (for the special history considered here and for any 
distribution p of random fields, using the Cauchy-Schwartz inequality) 

u 2l o(#i, H 2 f - u 2fi (Hi)u 2fi (H 2 ) < . (C35) 

In fact, for Hi ^ H 2 ( which excludes also Hi and H 2 being both equal to ±00, with the 
same sign), and R 7^ the expression in eq. (C35) is actually negative. This implies that 
for R 7^ the part in the quadratic action, which was coupling different magnetic fields, will 
be irrelevant at dimensions d < 6 under coarse graining. Its fixed point value is zero. In 
contrast, u 2 ,o(H) stays constant under coarse graining. 

What does this mean for the coarse grained action? If we leave out all the terms in the 
action that are zero or irrelevant at d = 6 — e, the different magnetic fields are completely 
decoupled, and the critical exponents (for R > R c at least) can be extracted from coarse 
graining the following action at fixed magnetic field H : 

S Ho = -J d d q J dtfj(-q,t)(-l/T d t + q 2 - X - 1 (H )/J)r ] (q,t) (C36) 
+ 1/6 J d d x J dtf)(x,t)(r](x,t) 3 u(H ) 
+ 1/2 J d d x J dti J dt 2 u 2 ^H ,H )ri(x,ti)fi(x,t 2 ). 

(The index and argument H serve as a reminder that all coefficients in this action are 
evaluated at the same magnetic field H . The time integrals extend from —00 to 00. The 
constant coefficients of the d t -tevm and the g 2 -term have been rescaled to 1 (see section 
VII A). 

APPENDIX D: BOREL RESUMMATION 

It is known [64,67] that the e-expansion yields only asymptotic series for the critical 
exponents as functions of dimension. It is important to note that asymptotic series in 
general do not define a unique underlying function. In the following we will give a definition 
of an asymptotic series, discuss several examples and give the special conditions under 
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which there does exist a unique underlying function. So far it has not been possible to 
show that the e-expansion would meet these conditions; we therefore cannot assume that it 
uniquely determines the critical exponents as functions of dimension. Applying the methods 
of reference [106] to our problem we Borel resum the perturbation series for our critical 
exponents to 0(e 5 ). A comparison of the result to the numerical values of the exponents in 
3, 4, and 5 dimensions is given in the main text. 

1. Definition of an asymptotic series 

Let the variable z range over a sector S in the complex plane with the origin as a limit 
point (which may or may not belong to S): 

Arg(z) < a/2, \z\ < \z \ . (Dl) 

A power series J2kLo fkZ k is said to be asymptotic to the numerical function f(z) defined on 
S as z — > 

oo 

f(z) ~ £ f k z k (D2) 

fc=0 

if the approximation afforded by the first few terms of the series is better the closer z is to 
its limiting value, which is zero in this case [107]. Formally: 

N 

\f(z)-J2f k z k \«\z\ N (D3) 

as z — > for every N (i.e. for any given e > there exists a neighborhood U t of the origin 
so that 

N 

\f(z)-J2f k z k \<e\z\ N (D4) 

k=0 

for all z common to U e and S). 

This is equivalent to the statement that the remainder of the sum after summing the 
first iV terms is of the order of the first neglected term and goes to zero rapidly as z — > 0: 
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(/(*) - E ~ f M z M (D5) 

V fc=0 / 

for every N, where /m is the first nonzero coefficient after f^. As shown in [107] this implies 
that there exists a constant C(N + 1) (i.e. a number independent of z) and a neighborhood 
[/ of the origin such that 

I (m - £ f k z k ) | < \C(N + I (D6) 

for all z common to U and S. 

Often the terms in the series at first decrease rapidly (the faster, the closer z approaches 
zero) but higher order terms start increasing again. If the series is truncated at the minimum 
term, one obtains usually the best possible truncated estimate for f(z) with a finite error 
E(z). If the coefficients C(k) are of the form 

C(Jfe) ~ C(k\fA^ (D7) 

with A > 0, (3 > 0, and C > real constants, the error E(z) can be estimated explicitly [67] 
to 

E(z) = min N C{k)\z\ k ~ exp[-/3(A/|z|)^] . (D8) 

which decreases rapidly as \z\ — > 0. It determines a limit to the accuracy beyond which it 
is impossible to penetrate by straightforward summation of a finite number of terms of the 
series. 



a. Non-uniqueness of f(z) 

One sees that an asymptotic series does not in general define a unique analytic function 
f(z) over S. Any other function g(z) which is analytic in S and smaller than E(z) in the 
whole sector can be added to f(z) such that in S the series Y^kLofk^ is asymptotic also 
to g(z) + f(z). If however the sector is big enough, i.e. if a > nj3 then according to a 
classical theorem (Phragmen-Lindelof) there exists no nonvanishing function over S that is 
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analytic in S and bounded by E(z) in the whole sector. In that case the underlying function 
is uniquely determined by the series [67]. 



b. Examples 

Let S A be the sector < \z\ < oo, = \Arg(z)\ < tt/2 - A with A > 0. 

1. The function exp(— 1/z) is asymptotic to Y^=q^ ■ z n — over Sa- 

|exp(-l/z) - 0| = exp(-cos(0)/|^|) < exp(- cos (tt/2 - A)/\z\) « \z\ N (D9) 
for any N as \z\ — > 0. 

2. Similarly one can show that over Sa the series 

oo 

(D10) 

is asymptotic to the functions 1/(1 + z^ 1 ^), (1 + exp(—l/z))/(l + z^ 1 ^) and 1/(1 + 



exp(— ^JJlJz)) + z^ ^) as z — > [107]. (A similar nonuniqueness is expected for the under- 
lying functions of the e-expansion.) 

3. The Stirling series is one of the oldest and most venerable of asymptotic series. It 
expresses the asymptotic behavior of the factorial function n\ for large values of n 



/27r\ 2 fn\ n . Ai A 



with Ai = 1/12, A 2 = 1/288, etc. It becomes the full asymptotic expansion of the Gamma 
function T(z) for complex argument z: 

i r w-(l) i (f)*( 1+ g^i«i(;) Vl -" (D12) 

for z — > oo and |^4r^(^) | < 7r [90]. 

2. Borel resummation 

The Borel summation is a technique for expressing an asymptotic series J2T=o fk zk as the 
limit of a convergent integral. For an introduction to this technique see reference [90,107]. 
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The idea is to modify the coefficients of the series such that for sufficiently small z it converges 
to a function B(z), which can be calculated. In the simplest case 

oo 

B{z) = Y.h/k\z k (D13) 

k=0 

would be appropriate. Suppose that also 

f(x) = / exp(-t)B(xt)dt (D14) 
Jo 

exists. Then, by expanding the integral / °° exp(—t/x)B(t)dt/x for small x and integrating 
term by term, one finds [90] that f(x) = DfcLo fk xk as a; + . 

a. Example 

The series J2 fkZ k = J2^=o n\{—x) n diverges, but B(x) = Y^=o(~ x ) n converges for \x\ < 1 
to 1/(1 + x). Thus the Borel sum of E^n^-a;)" is f( x ) = Jo°° <*p(-*)/(1 + xt ) dt - For 
z in the complex plane, f(z) is asymptotic to XX— z) n n\ over the sector S: < \z\ < oo, 
| At g(z) | < Ti — e, e > 0. Since \f n \ ~ n\ we have j3 = 1 from eq. (D7), and since a > n 
for small enough e, the uniqueness theorem of section D 1 a implies that f(z) is the unique 
underlying analytic function of the series over S. 

3. Borel resumming the e-expansion 

It is known [67], that the e-expansion series is an asymptotic series. It has a zero radius 
of convergence. Lipatov, Brezin, LeGouillou, and Zinn- Justin have shown [108], that the 
coefficients of higher orders of the e expansion /(e) = J2k(~ e ) h fk (with e = 4 — d) have the 
asymptotic form fk ~fe^oo k\a k k b c. (The factorial growth of the coefficients reflects the zero 
radius of convergence.) Here, / stands for r\ or \/u (or uj or g , see [106]), and a = 1/3, 
b = 7/2 for 77, and a = 1/3, b = 9/2 for l/v. 

Using the newest results for the coefficients to 0(e 5 ) as given in Kleinert et al. [59] and 
the asymptotic form given above, we resum the coefficients to 0(e 5 ), using the method 
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discussed in reference [106,109] and reference [110]. The method is based on a modified 
Borel transformation: 

/>00 

/(e) =/ dx/(ea)exp(-x/(ea))(x/(ea)f +3/2) B(x). (D15) 

JO 

One obtains B(x) = Ek(-x) k B k , with B k = f k /(a k T(k + b + 5/2)) c/k 3 / 2 . 

B(x) is given by its Taylor series only for \x\ < 1. It can be continued analytically beyond 
the unity convergence circle in the complex plane, using the conformal mapping [111] x — > w 
with w(x) = [(1 + x) 1/2 - 1]/[(1 + x) 1/2 + 1], or x = Aw/(1 - w) 2 . The integration interval 
[0, oo) then goes over into the interval [0, 1] and the cut (— oo, —1] goes over in the boundary 
of the unit disk. The expansion of B(x(w)) in terms of w converges to B(x) in the entire 
region of integration w G [0,1). The coefficient of w N is determined on the basis of the 
coefficients f k for k < N of the initial e expansion. Neglecting the terms of 0(w m ), with 
m > N, one obtains an approximation //v(e) which takes into account the 0(e N ) corrections 
to /. 

Since we are working with truncated rather than infinite series there is an arbitrariness 
in the reexpansion of B(x) in terms of w [110,106]. Rather than expanding B(x(w)) di- 
rectly, Vladimirov et al. actually expand ((1 — w) 2 ) x B(x(w))/A = B x (w) in terms of w, 
so that (since x = Aw /(I — w) 2 ) the result for B(x(w)) is given by B(x) — > x x B x (w) = 
(x/w) x Ek B[ k) w k . X is a parameter that must be chosen correctly, such that /jv(e) is 
matched to the asymptotic form /(e) as e — > oo, given by: /(e) ~ e -_ >00 (e) A . Since the 
asymptotic form of the exponents relative to e is unknown, it is argued in [109] and [110] 
that A should be fixed from the condition of fastest convergence of /at, i.e. it should mini- 
mize A at = |1 — /jv(1)///v-i(1)|. We have written a computer program that performs this 
procedure for the exponents \jv and r\ using the following results for the e-expansion to 
0(e 5 ) from reference [59]: 

\fv = 2 - e/3 - 0.1173e 2 + 0.1245e 3 - 0.307e 4 + 0.951e 5 + 0(e 6 ) , (D16) 
r] = 0.0185185e 2 + 0.01869e 3 - 0.00832876e 4 + 0.02566e 5 + 0(e 6 ) . (D17) 
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Through the perturbative mapping of our problem to the pure Ising model in two lower 
dimensions to all orders in e these are also the results for the exponents in our system in 6 — e 
dimensions. Figure 14 shows the comparison between the Borel resummed exponents r], 
and (38 — |(d — 2r] + fj), and the corresponding numerical results at integer dimensions. We 
have used that in perturbation theory fj — rj, as follows from the mapping to the equilibrium 
random field Ising model [60,62]. 

The standard Borel resummation has been modified in various ways specifically to resum 
the 4 — e expansion for the exponents of the pure Ising model. Le Gouillou and Zinn- Justin 
[110] have developed a resummation prescription that takes into account that the critical 
exponents have a singularity at the lower critical dimension d = 1 (i.e. e = 3). In a 
more recent calculation [110], they treated the dimension of the singularity as a variational 
parameter to improve the apparent convergence of the expansion. They also imposed the 
exactly known values of the exponents in two dimensions. Their results agree (partly by 
design) very well with the numerical exponents of the pure equilibrium Ising model in two 
and three dimensions. Below two dimensions the apparent errors of their results are rapidly 
increasing. 21 

In our problem the lower critical dimension is probably 2 rather than 3, i.e. e = 4 rather 
than e = 3. In our Borel-resummation of the 6 — e expansion results for our critical exponents 
we have made no assumptions about a singularity for the critical exponents, nor did we 
impose any other information. Instead, we have applied the conventional Borel-resummation 
as adapted by Vladimirov et al, which shows good agreement with our numerical data (see 
figure 14). Our resummation of corrections for v up to order 0(e 5 ) has a pole around 2.3 
dimensions. (In reality v is expected to diverge at the lower critical dimension.) The pole 
did slowly increase with the order to which the resummation was performed. It is not clear 



We should note that their work was based on a form for the fifth order term which later turned 
out to be incorrect [59]. 
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however whether this tendency persists at higher orders. 

Figure 19 shows a comparison of our Borel resummation for v with the results of Le 
Guillou and Zinn- Justin from 1987 [110]. It is believed that the "straightforward" Borel- 
resummation according to Vladimirov et al. [106,109] and the resummation according to 
LeGuillou and Zinn- Justin [110] which treats the dimension of the singularity in the expo- 
nents as a variational parameter, should both at high enough orders converge to the pure 
Ising exponents, although this has not been proven. Unfortunately, the curve from LeGuillou 
and Zinn- Justin [110] shown in figure 19 has been obtained using an epsilon expansion for 
the critical exponents, which later turned out to be incorrect in the 5th order term [59]. Our 
own Borel-resummation however, which is also shown in figure 19, is based on the newest, 
presumably correct result for the fifth order term [59]. We are currently trying to duplicate 
the variable-pole analysis of LeGuillou and Zinn- Justin with the correct fifth order term. 

As we have already mentioned, there is no known reason to assume that the e-expansion 
is Borel-summable, i.e. that the values of the exponents would be uniquely determined 
as functions of dimension by the Borel summation of their asymptotic e-expansion. At 
this stage we do not know the size of the section of analyticity for the various Borel sums 
[110]. Hence we expect the e-expansion to give the correct value of the exponents only 
asymptotically as e — > 0, even after Borel resummation. There is no reason to expect that 
its extrapolation to e = 3 can be used to determine the lower critical dimension, i.e. the 
dimension where v has a pole, and where the transition disappears or loses its universality 
(see also section VIII B in the text). 

APPENDIX E: INFINITE AVALANCHE LINE 

In most of this paper we have focussed on the critical endpoint (R c , H C (R C )), in particular 
as it is approached from R > R c . In this appendix we discuss the onset of the infinite 
avalanche for R < R c . Our e-expansion can be applied to the entire line H C (R), R < R c at 
which the infinite avalanche occurs (with some reservations which we will discuss later). In 
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mean-field theory the approach to this line is continuous with a power law divergence of the 
susceptibility x ~ dM/dH and precursor avalanches on all scales (see appendix A). From 

Sh = - £ / dtJ^Jfj^tMt) - Y.J dtfuW-dt/Fo - <*] Vj (t) (El) 
+£(1/6) / dt%(*)Mf)) 3 « + (l/2)£ / *i / ^2%(ti)%(t 2 Ko, 
quoted from in eq. (65), and from the rescaling of the vertices given in eq. (82) 

<n = fe[-( m+ ") +2 ^ 2+2 " Mm ,„ . (E2) 

one sees that on long length scales the effective action is purely quadratic above 8 dimensions. 
This suggests that there is a continuous transition (as H approaches H C {R)) with mean- 
field critical exponents and a diverging correlation length £(x _1 ) with the scaling behavior 
^(X- 1 ) = b^x' 1 ), i.e. i ~ (X- 1 )- 1 ' 2 - Since X ~ X ~ - H C {R)\ (see appendix A) it 
follows that £~\H- H C (R)\- Uh with i/ h = 1/4 for d > 8. 22 

For d = 8 — e (e > 0) the vertex u> in the action Sh becomes relevant. In contrast to 
the critical endpoint where x' 1 — and w — 0, the infinite avalanche line is characterized 
by the "bare values " x~ l = and w = -2J 2 p'(-JM(H c (R)) - H C (R) + k) ^ 0. Figure 20 
shows the lowest order the correction to vertex w. With the Feynman-rules of appendix C 2 
and section VII B 1 the recursion relation to the same order becomes 



w'/2 = b^ 2 ^ |«,/2 + Ko/2) (W2) 3 8/(47r) 4 £ <kVtf - X' 1 1 Jf\ ■ (E3) 



Performing the integral over the momentum shell A/b < q < A and writing M d / 2 + 4 ) 
6(«/ 2 ) = l + e/2 log 6 we find 



22 If the Harris-criterion is not violated through the presence of large rare nonperturbative fluctu- 
ations in an infinite system, such as a preexisting interface (for a discussion see appendix G), i.e. 
if Vh > 2/d is a valid exponent-inequality, then the mean-field critical exponents with = 1/4 are 
only correct for d > 8, which is consistent with our result from perturbation theory. 
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w'/2 = w/2 + (w/2) (e/2 + u 2fi (w/2) 2 4/(47r) 4 log ft) . (E4) 

Since -u 2 ,o > 0, this equation has only two fixed points w* with w' = w for e > 0: either 
u>* = or w* = oo. Any system with bare value w^O will have effectively larger w on longer 
length scales. The system flows to the strong coupling limit. We interpret this as indication 
that in perturbation theory the transition is of first-order type below 8 dimensions. Indeed, 
our numerical simulation suggests an abrupt onset of the infinite avalanche in 2, 3, 4, 5 and 
7 dimensions, but smooth in 9 dimensions [45]. 

As we discuss in appendix G there are some questions as to whether in an infinite system 
the onset of the infinite avalanche would be abrupt in any finite dimension due to large rare 
preexisting clusters of flipped spins which provide a preexisting interface that might be able 
to advance before the perturbatively calculated critical field H C (R) is reached. These large 
rare fluctuations might be nonperturbative contributions which are not taken into account 
by our e-expansion. The progression of a preexisting interface has been studied previously 
in the framework of depinning transitions [30,29,28]. Our preliminary numerical simulation 
results in 9 dimensions do however seem to show a continuous transition at the onset of the 
infinite avalanche as predicted by the RG calculation [45]. 

APPENDIX F: DETAILS FOR THE e-EXPANSION OF THE AVALANCHE 

EXPONENTS 

1. The second moment of the avalanche size distribution 

In this section we show that the second moment (S 2 ) / of the avalanche size distribution 
D(S,r, h) scales in the adiabatic limit as 

(S 2 ) / ~ / dti / dt a dt p d d x a d d xp(s a (t a , x )s a (t , x a )s p {t p , x )s l3 (t 1 ,x (3 )) f . (Fl) 
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where a and (3 specify the corresponding replica that have identical configurations of random 
fields and are exposed to the same external magnetic field H(t) — H + £lt, with f2 — > 0. 23 
We start by computing the (not yet random field averaged) expression 

J dt a d d x a dtpd d x t3 {s a {t a) XQ)s a {t 0) x a )s p {t ()) XQ)s p {t 1 ,x (5 )} (F2) 
= J dt a d d x a {s a (t a ,x )s a (t ,x a )} J dtpd d x f3 {s l3 (tp,x )s l3 (t 1 ,xp)} 

where { } stands for the path integral over the product with the ^-function weight in Z that 
singles out the correct path through the space of possible states for the given configuration 
of random fields and the given history In eq. (F2) the two replicas are uncoupled since we 
have not yet averaged over the random fields. As we have seen in appendix C 

(AS/AH) a = J dt a d d x a {s a (t a ,x )s a (t ,x a )} = (F3) 

/o 
dt a d d x a — lim AH ^ 
Ota 

(s a (t O) a:a)|H« o (to)=H(t o )+A J f/0(to-t C( ) - s(to,X a )\ HSo (to)=H(to)) / AH 

is the response of replica a to a perturbing pulse of amplitude AH applied at field H(t a ) at 
site x integrated over the entire system. 

If no spin flips in response to the perturbation, the total response will be 

(AS/ AH) = AS harmomc /AH = C 2 (F4) 

where C 2 is a constant that depends only on the parameters k, J, and the coordination 
number z of the lattice. 

If, on the other hand the perturbation triggers an avalanche of spin flips from the lower 
to the upper potential well, AS = Sfl ip = S a will be of the order of the number of spins 
participating in the avalanche (see also appendix CI). 24 



23 A heuristic justification of this was given in section IX together with an explanation of why 
replicas are necessary. 

24 S a is not exactly equal to the number of spins flipping in the avalanche. It contains also the 

98 



The expression in eq. (F2) is the product of the total response to the same perturbation 
at site x measured in replica a at time t and in replica (5 at time t\. At finite sweeping rate 
fl/k the corresponding values ((AS) a / AH) and ((A^^/Aif ) do not have to be the same, 
since the responses are measured at potentially different values of the external magnetic field 
(H = H(t ) and Hi = H(ti) respectively). (We only consider the adiabatic case, in which 
the sweeping rate fl/k is small compared to the relaxation rate T k, so that the magnetic 
field can be assumed to be constant during the course of an avalanche. 25 ) 

Without loss of generality let us assume that H(t±) > H(t ). First we discuss the case 
that there is an avalanche S a triggered by the perturbation of amplitude AH in replica a at 
field H . We further assume that t\ is much bigger than t , such that Hi > H + AH. In this 
case the response to the pulse in replica (5 will be substantially different from the response 
S a in replica a. The spins that are pushed over the brink by the perturbation at field H 
in replica a, will in replica (3 be triggered by the increased external magnetic field before it 
reaches the bigger value Hi at which the response is measured. For fl/T , AH, and Hi — H 
small enough, the response in replica (5 at field Hi will then be just the harmonic response 



harmonic response that each spin flip causes through the coupling to the neighboring spins. This 
harmonic response couples back to the original spin and propagates to the next-nearest neighbors 
with an amplitude damped by the factor Jij/k and so on. Occasionally it may cause an avalanche 
to continue which would otherwise (in the hard spin model) have come to a halt. However since 
this is a short-ranged effect, we do not expect it to be of any relevance to the scaling behavior on 
long length scales. In mean-field theory the harmonic response only amounts to a constant factor 
relating S a to the number of spins participating in the avalanche. 

25 We take the adiabatic limit fl — ► at finite correlation length £, before approaching the critical 
point of diverging avalanche size and time to avoid triggering a new avalanche before the previous 
one has come to a halt. This is consistent with our computer simulations at finite system sizes 
where avalanches occur only sequentially. 
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C 2 or a different avalanche. If it is the harmonic response, the expression in eq. (F2) takes 
the form 

(AS/AH) a (AS/AH)p = {SJAH)C 2 . (F5) 

Similarly one might imagine scenarios in which there is an avalanche Sp triggered only in 
replica (5, it i.e. 

(AS/AH) a (AS/AH)p = (Sp/AH)C 2 , (F6) 
or where there is no avalanche happening at either field value 

(AS/AH) a (AS/AH)p = (C 2 ) 2 . (F7) 
It is also possible that two different avalanches are triggered in the two replicas: 

(AS/AH) a (AS/AH)p = (S a /AH)(Sp/AH) (F8) 

with S a ^ Sp. 

We are interested however in contributions due to the same avalanche response in both 
replicas 

(AS/AH) a (AS/AH)p = (S a /AH)(Sp/AH) (F9) 

with S a = Sp. As we have seen, a necessary condition is that H — AH < Hi < H +AH. We 
denote with Pfl ip = c AH + o((AH) 2 ) (with c a constant in the critical regime) the fraction 
of all possible configurations of random fields in which a local perturbation of amplitude AH 
at field H, applied at site x , causes at least one spin to flip. For Q and AH small enough 
the fraction of all possible configurations of random fields in which the local perturbation 
will lead to the same initial spin flip triggering the same avalanche S in replica a and replica 
f3, is to leading order in AH proportional to the size of the overlap P^ 1 of the two intervals 
[H , H + AH] and [H u H x + AH], multiplied by P flip , with 

PflT = (1 - e(l#i - H \ - AH))(AH - \H ± - H \)/AH (F10) 
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(see figure 21). We can now compute the random-field average of the expression in eq. (F2), 
denoted by ( )/ to leading order in AH 

AS AS \ _ r< (S 2 ) f both , r, ( s )f (t P both\ p 



+ (C 2 ) 2 + (S a S,)/(AH) 2 P 2 fltp (l - P^) , (Fll) 



where {S 2 )f is the mean square avalanche size, and (S)f is the mean avalanche size, and C\ 
and C2 are constants in the critical regime. The last term accounts for cases in which two 
different avalanches S a 7^ Sp are triggered in the two replicas. 

The last three terms in eq. (Fll) approach a constant as AH — > 0, since Pfn p ~ AH. We 
will now analyze the first term, which is proportional to (S 2 ) in more detail. The function 
multiplying (S 2 ) j is sharply peaked around H = Hi (see figure 21). Since Pfl ip ~ AH it is 
proportional to P^ h / AH . From eq. (F10) we have 

/ dffiPjg* AH = 1 (F12) 

J H -AH 1 

independent of AH. The same integral applied to the other terms in eq. (Fll) yields 
contributions of order O(AH) which are negligible compared to the first term as AH is 
chosen small. With H 1 = H + £lt we can express the integral in terms of time 



r AH/n 

/ QdhP^/AH = 1 . (F13) 
J-AH/n J 



We then obtain 



rAH/n AS AS 

lim AH ^ lim n ^ / £Mi( A „ ) f = C^S 2 ) f . (F14) 

J-AH/n Ail a 

With eq. (F3) this leads to the scaling relation 

(S 2 ) f ~ JdhJ dt a dtpd d x a d d xp(s a (t a , x )s a (t , x a )s p {t p , x )sP(t u x p )) f . (F15) 

which was to be shown. In this notation we have suppressed the factor Q and the various 
limits for clarity. The integrals over time extend from —00 to +00 with an infinitesimal 
associated change in magnetic field. 
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2. Feynman rules for two replicas 



We study the behavior of S a(3 of eq. (161) under coarse graining analogously to the 
calculation done before for just one replica, with the difference that instead of two, there 
are now four fields to be considered (two for each replica). 26 In the following section we 
briefly describe the associated Feynman rules. This section may be skipped by the reader 
uninterested in the details, since it turns out that there are no loop corrections to 0(e) to 
(S 2 )- 

In the Feynman graphs for the loop corrections, the fields of the a replica are symbolized 
by arrows on full lines, whereas those for the (3 replica are symbolized by arrows on dashed 
lines. A vertex u mnpq has then m outgoing arrows on full lines, n incoming arrows on full 
lines, p outgoing arrows on dashed lines and q incoming arrows on dashed lines. In this 
notation, the fact that u 0npq = if n ^ 0, u mn0q = if q ^ and u 0n0q = 0, which we 
discussed in section IX B, means that any vertex with incoming arrows of a certain replica 
must have at least one outgoing arrow of the same type of replica, i.e. there are no "sinks", 
with only incoming lines of a certain replica. Furthermore, since the spins from different 
replica do not interact directly, and since wo,i,i,o — M i,o,o,i = 0, there are only two kinds of 
propagators, one for each replica. In other words, in any diagram, an outgoing line can be 
connected only to an incoming line of the same replica. 

Using the above rules and causality, one finds that corrections to vertices with lines of 
only one replica, can only receive corrections from vertices of the same replica. There are 
no contributions from diagrams that also involve the other replicas. That means that our 
results for the magnetization and other quantities that can be calculated using only one 
replica, are unaffected by the introduction of a second replica. 

"Pure" (or one-replica) vertices which depend on more than one time usually have several 



26 In section IX B we already derived the appropriate partition function. In this appendix we use 
the notation introduced there. 
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different contributions. For example the vertex 1*2,2(^1, fa, £3, £4) has two main contributions 
that are obtained by partial integration of the corresponding term in the action as discussed 
in appendix C. One contribution is derived from 1x2,2(^1, fa, fa, fa) and has fa = and t 2 = t%. 
The other contribution is derived from U2,2(fa, fa, fa, fa) and has fa = t£ and t\ = t%. If we 
have two replicas, there are different "mixed" vertices each of which corresponds to one such 
the individual contribution to a pure vertex. Each of them has the same bare value as such 
a pure counterpart, since both are obtained in the same way from mean-field theory. With 
"corresponding" we mean that the times associated with the different legs of the mixed 
vertex, are assigned in the same way to the legs of the corresponding part of a pure vertex. 
The part of a pure vertex formally corresponding to 1x1,1,1,1, for example, is given by that 
contribution to 1x2,2(11, fa, fa, fa), which has fa = t~£ and fa = t%. Conversely the part of a pure 
vertex corresponding to 1*1,2,1,0 is given by that contribution to 1x2,2(^1, fa, fa, fa), which has 
fa = £3 and t\ = £4 . Notice, that in any mixed vertex all legs carrying a certain time label 
(one outgoing and any number of incoming arrows), must belong to the same replica. In 
contrast to the different contributions to one pure vertex, the corresponding mixed vertices 
do not add up to a single mixed vertex, since they are multiplied by field from different 
replicas rf 7^ rf . 

The loop corrections to mixed vertices formally look the same as those the corresponding 
parts of the pure vertices. For each loop correction to a mixed vertex there is a matching 
correction to the corresponding part of the pure vertex and vice versa. The combinatoric 
factors are also the same. This implies in particular that choosing the same spin rescaling 
for both replica as we did before in the case of only one replica, renders marginal not only 
1*2,0 and 1X2,0, but also 1*1,0,1,0- 

3. Scaling of the second moment of the avalanche size distribution 

We need to find the scaling behavior of the "Greens function" 

(s a (t a , x )s a (t , x a )s p (t^, x )s l3 (fa, Xp)) f 
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(F16) 



from its behavior under coarse graining. The topology of the diagrams permits no 0(e) loop 
corrections to the corresponding vertex function. 

One finds the canonical dimensions of the fields [105] (where "dimension of" is denoted 
by "[ ]" and A is the upper cutoff in momentum): [r](p,iu)] = A~ d / 2 ~ 2-2 , [f)(p,uj)] ~ A~ d / 2 . 

For calculating Greens functions one introduces source terms in the action. From the 
(functional) derivative with respect to the source fields, one obtains the corresponding aver- 
age correlation functions. In the end the source fields are taken to zero again, since usually 
they have no physical significance. In our case the following three source terms are needed: 
/ d d q J duL(q, ou)r)(q, u>), J d d q J duL(q, u>)f)(q, cu), and the term needed for the calculation of 
the (spacially) composite operator in {S 2 )f, given by 

J d d q J dui j dw 2 L 2 (q,tu u u 2 ) J d d qf)(q,Ui)f)(p - q,u 2 ) ■ (F17) 

L, L and L 2 are the respective source fields: the corresponding canonical dimensions are 
[L(q,u)] ~ A' d / 2 + 2 , and [5/5L(q,u)} ~ \{d/2-2) K (-d-z) _ A -d/2-2-,_ Similarly [L(q,u)] ~ 
A-<*/2-^ and [S/SL{q,uj)] ~ A^ 2 . And also [L 2 (p, O)] ~ A^ d ~ 2z \ and [5/5L 2 {p,Q)] ~ A . 
From eq. (164) and the fact that Greens functions in the fields r\ and fj scale in the same way 
as those in terms of s and s (see section VI B), we then find (without loop corrections) that 
{S 2 )f ~ A~( 4+z \ Below the upper critical dimension, the canonical dimensions of the fields 
r](q,u) and fj(q,u) are corrected by A^/ 2 ^ and A ( ^~ f ' //2 ) respectively. With r] = fj from the 
mapping to the pure Ising model [60], one obtains (to 0(e)) (S 2 ) / ~ A _ ^ + ( 2_ ^ 2 ^. Similarly, 
one finds for the higher moments (S n ) f ~ ^-((n-i)z+(2-v)n) to o(e). In section IX D this 
result is compared to the scaling behavior of (S 2 ) as obtained from the scaling form of the 
avalanche size distribution 

{S 2 } rvrW'SPih/rP*) (F18) 

(with the appropriate scaling function «S±) to extract the results for 1/sigma and r. There 
we obtain 

l/a = 2 + e/3 + 0(e 2 ). (F19) 
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and 



r = 3/2 + 0(e 2 ). 



(F20) 



APPENDIX G: RELATED PROBLEMS 
1. Comparison with depinning transitions 

a. Relation to fluid invasion in porous media 

Much progress has been made in the study of fluid invasion of porous media [112,27,28]. 
A preexisting fluid in the porous medium, for example oil in rock, is displaced by an invading 
fluid (water), which is driven by an applied pressure P. The interface between the two fluids 
is pinned at pressures lower than a threshold value P c and advances continuously at higher 
pressures. The interface depinning transition is accompanied by a diverging correlation 
length, critical fluctuations and universal exponents. There are three different universality 
classes for the associated critical exponents, corresponding to low, intermediate and high 
disorder in the system. At low disorder the marginally stable interface at P c is faceted, at 
intermediate disorder it is self-affine and at high disorder it is self-similar. 

Already before we did our work on hysteresis, Robbins, Cieplak, Ji and Koiller have 
pointed out the analogy of fluid invasion to the physics of domain wall motion in Ising 
ferromagnets [28,30,29]. There, the interface separates regions of up and down spins. As 
the magnetic field is ramped, one domain grows at the expense of the other — the interface 
is pushed forward. Quenched disorder may be due to random fields or random bonds. The 
authors study the zero temperature nonequilibrium RFIM with a rectangular distribution 
of random fields of width A and a preexisting flat interface separating up spins from down 
spins. They increase the external magnetic field adiabatically and study avalanches of spin 
flips that are triggered by flipping spins at the interface. As in our system a spin in such 
an avalanche flips up if its local effective field becomes positive. The main difference is that 
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in their system no spins are allowed to flip ahead of the interface. Each spin flip can be 
interpreted as an advancement of the pre-existing system-spanning interface. The value of 
the magnetic field at which any part of the interface reaches the other side of the system is 
the critical field which corresponds to the threshold pressure P c in fluid invasion in porous 
media. As in fluid invasion, Ji and Robbins find in 3 dimensions a faceted regime at low 
disorder A, a self-affine regime at intermediate A and a percolation regime at large A. At 
a critical width A2 = 3.41J, which separates the self-affine regime from the self-similar 
regime, a diverging length scale in a bulk quantity ( "fingerwidth" ) is observed. The critical 
exponents associated with the interface depinning transition in the self affine regime have 
also been calculated in an e-expansion [30,29]. They are the same for random-field and 
random-bond disorder [29]. 

We expect that in our system it would be possible to observe the same critical behavior 
at the onset of the infinite avalanche in large enough systems with less than critical disorder 
(R < R c ). If the system is big enough, there will certainly be somewhere a rare large cluster 
of flipped spins, even at relatively low magnetic fields. As the field is slowly raised, the 
surface of such a cluster is expected to act as a preexisting interface analogously to Ji and 
Robbins system. The small clusters that are flipped ahead of the interface probably have 
no influence on the long length scale behavior and the critical exponents associated with 
interface progression such as roughening exponents. The onset field for the infinite avalanche 
in an infinite system corresponds to the threshold field at which a preexisting interface gets 
depinned in Ji and Robbins system. 

Numerical and analytical studies [19,20,29,30,28,21,22,31,23] leave no doubt that the 
interface depinning transition has an associated diverging height-height correlation length 
and critical fluctuations. There are scaling forms for quantities related to the shape of the 
interface and its progression, which suggest that the interface gets depinned in a second 
order transition. Nevertheless in our system (in three dimensions) we have called the onset 
of the infinite avalanche in systems with less than critical disorder an abrupt or 11 first order" 
transition. What's come over us ? Our numerical simulations in three dimensions clearly 
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show a kink in the magnetization curve at the threshold field H C (R), rather than power law 
behavior as expected in continuous transitions. This is not a contradiction. One has to be 
careful about which quantities are measured: The critical fluctuations at the interface yield 
no contribution to bulk quantities (e.g. the magnetization), which are measured per unit 
volume. In references [126,128] we propose an experimental setup that would allow one to 
observe bulk and interface fluctuations simultaneously. In the conventional setup for mea- 
suring magnetic hysteresis loops discussed in this paper however, the second order depinning 
transition with a diverging height-height correlation function at the interface discussed by 
Robbins et al. is buried invisibly inside the line of abrupt ("first order") transitions seen 
in our magnetization curves in sufficiently low dimensions. Note that this connection can 
only be established for large enough systems, where rare large fluctuations provide a preex- 
isting interface of flipped spins. In smaller systems, disorder induced nucleation effects will 
determine the size of the onset field for the infinite avalanche in our hysteresis model. 

At the critical disorder Ag, separating the self-affine regime from the self-similar regime, 
Ji and Robbins find a diverging bulk length scale, the "fingerwidth" of regions of flipped 
spins. The critical exponents associated with this transition are different from the exponents 
seen at our critical endpoint (R c , H C (R C )). This is not too surprising: unlike the case for 
R < R c , near R c a large advancing interface in our problem runs into pre-existing flipped 
regions of all sizes — presumably changing the universality class. Also we have several 
infinite fronts at R c . Ji and Robbins report that for a bounded nonanalytic distribution of 
random fields the corresponding exponents are nonuniversal in the depinning problem, but 
depend on the analytic form of the edges of the distribution of random fields. However, 
for an analytic distribution they presumably take universal values [113]. In our problem we 
use an analytic distribution and v is universal, as is shown by the RG calculation. It has 
been shown [57] however that the exponents are different for a rectangular distribution of 
random fields in our system also. Perhaps it would not be universal in our problem either 
for rectangular distributions of random fields [57] (or other distributions with singularities 
at their tails). In section G2 this issue will be discussed in more detail. 
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Above a critical amount of disorder, using single interface dynamics, Ji and Robbins see a 
self-similar regime, where the interface grows to a percolating cluster at a certain "threshold 
field" H^ m9le (R). In our dynamics in an infinite system such a percolating cluster might 
occur already at a lower value of the external magnetic field, because it is probably easier 
to connect preexisting clusters of flipped spins to form a percolating cluster than it is to 
push a domain wall through the system all the way till it reaches the the opposite side. 
One would expect the same argument to apply for R < R c , i.e. H C (R) < H^ ngle (R). In 
references [126,129] we discuss four potential experimental setups for measuring Barkhausen 
noise (or avalanche size distributions in general), explaining for each case whether the single 
interface model of Ji and Robbins or our model with many interfaces and domain nucleation 
is expected to apply. 

2. Other models 

a. Adiabatic models for hysteresis 

There are several numerical studies of related hysteresis models, which we discuss in 
more detail in reference [126,128]. Here we only list a few closely related models. 

Maslov and Olami [57] have simulated the same model as we have studied, but for a 
rectangular distribution of random fields [57] rather than a Gaussian. The authors find 
different critical exponents than ours in 3 dimensions and in mean-field theory. They report 
that their numerical results speak in favor of the upper critical dimensions being d c = 6. 
They claim that the model does not belong to any known universality class. 27 



27 Presumably in each dimension there will be a critical power law for the tails of a bounded 
distribution, so that distributions with a power law larger than this critical value will lead to 
different critical exponents. In this sense the tails of a rectangular distribution can be thought of 
as an infinite power law, and conceivably might lead to a different universality class than ours in 
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Koiller, Ji and Robbins [28] point out that for a rectangular distribution at R < R c there 
is a close connection to diffusion percolation and bootstrap percolation in 2 and also in 3 
dimensions [113], which are modified percolation models in which the occupation of a site 
depends on its environment [114]. 

Coram, Jacobs, and Heinig [115] have studied the zero temperature nonequilibrium ran- 
dom bond Ising models, both spin glasses (SG) and random ferromagnets (RFM) with 
nearest neighbor interactions in one, two and three dimensions. They report power law 
sensitivity to single-spin-flip perturbations in 2 and 3 dimensional random ferromagnets, if, 
starting from the spin state with all spins pointing down (ground state at H = 0) the field 
has been raised to a positive critical value. The associated critical exponents are r = 1.37 
in 2 dimensions, and r = 2.8 in three dimensions. It is likely that the value of the exponent 
r would be different for a Gaussian distribution of random bonds (see footnote 27 on page 
108). 

Vives et al. [50] have studied that case for a negative mean of the Gaussian bond dis- 
tribution. They trigger avalanches by a slowly increasing homogeneous external magnetic 
field (as in our model). For the avalanche size distribution integrated over one branch of the 
hysteresis curve they find r + a (35 = 2.0 ± 0.2 in three dimensions and r + a f38 = 1.45 ±0.1 
in two dimensions, which are rather close to the values of the nonequilibrium RFIM (see 
section G2c). One would expect the value the exponent r to be smaller than r + a (35, i.e. 
their result seems to deviate from the number obtained by Coram et al. for a rectangular 
distribution of ferromagnetic bonds with strengths between and 1. That is not surprising 
if we consider the differences in the two approaches. It might be interesting to compare the 
shape of the hysteresis loops of the two models at various values of the relevant tunable 
parameters. A similar comparison of our hysteresis loops to those of Maslov and Olami's 
RFIM with a rectangular distribution of random fields [57] revealed marked differences of 



all dimensions. 
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the two models. For example, in our model at less than critical randomness, there were 
precursors to the infinite avalanche, while in Maslov and Olami's model there were none. 

Bertotti and Pasquale [118] have studied hysteresis phenomena in the Sherrington- 
Kirkpatrick spin glass model [117] with N Ising spins Sj = ±1 on a lattice with random 
infinite range interactions that are distributed according to a Gaussian with mean zero. 
They used a slightly generalized dynamics that allows for temperature like fluctuations at 
a finite sweeping frequency Q with H = Qt. They report power law scaling of the power 
spectrum due to avalanches of spin flips near the central part of the saturation hysteresis 
loop, where the rate of change of the magnetization dM/dt is approximately constant. The 
authors note that this scaling behavior as well as the associated shape F(uo)/VL very much 
resembles the behavior observed in Barkhausen effect experiments in soft magnetic materials 
[119]. 

For sweeping rate Vt — > we can think of their model without random fields as an 
infinite range mean-field theory for the nonequilibrium random bond Ising model, simulated 
by Vives et al. [50]. Vives et al. found that in two and three dimensions the nonequilibrium 
RFIM does reveal a critical point of the kind we found in the nonequilibrium RFIM. It would 
be interesting to look for the same kind of critical point in Bertotti et aVs nonequilibrium 
SK- model (with and without random fields), for example by setting the mean J of the 
distribution of random bonds to a nonzero value and tuning the widths of distributions of 
random fields and random bonds. 

Rudyak [37] has suggested a theory for dielectric hysteresis in ferro electrics which leads 
to the same picture as our hard-spin mean-field theory for R < R c . 

b. Dynamical and other hysteresis models 

There are several studies of scaling behavior of the area of the hysteresis loop as a 
function of driving frequency Q and amplitude H Q of the of the external magnetic field 
H(t) = if sinf2t, and other dynamical effects [120,118,121]. Most of these model contin- 
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uum and lattice spin systems do incorporate temperature effects, but no quenched disorder. 
Rao, Krishnamurthy and Pandit give a nice review and discussion of previous (empirical) 
hysteresis models such as the Preisach model and mean-field types of theories like the Stoner- 
Wohlfarth theory and others [121]. 

c. Conjectures about other models in the same universality class 

Recently Vives, Goicoechea, Ortin, and Planes found [50] that the numerical exponents 
v, j3, t, and z in the nonequilibrium zero temperature RFIM and the RBIM (random bond 
Ising model) have very similar values in two and three dimensions. In two dimensions, the 
exponents for the RFBEG (random field Blume- Emery- Griffiths model [122]) seem to be 
similar also. In this interesting paper the authors suggest that these models might actually 
be in the same universality class. Admiring their work, we however have some concerns as 
to whether their critical exponents will remain unchanged for larger system sizes: they used 
systems of linear size up to L = 100; we used much larger systems, up to 7000 2 and 800 3 for 
the RFIM, and found that finite size effects are actually quite prominent and lead to shifted 
results for the exponents. Although the equilibrium versions of these models are not in the 
same universality class, it is known [30,29], however, that the nonequilibrium single interface 
depinning transitions of the RFIM and the RBIM do have the same critical exponents. 

In the following section we will discuss some symmetry arguments, that would indeed 
speak in favor of Vives et aVs conjecture and would suggest that the universality class 
of our model extends even beyond just the RBIM. A large universality class would also 
explain the surprisingly good agreement with experiments discussed in section X. Generally 
one may ask how robust the universality class of our model is against the introduction of 
other kinds of disorder, other symmetries for the order parameter, long range interactions, 
different dimensions, and altered dynamics. The variation with dimension has already been 
discussed at the appropriate places in this paper (see for example section X). 
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Other forms of disorder and symmetries 

If a new kind of disorder in an otherwise unaltered system changes neither the symmetries, 
nor the interaction range, nor the dynamics, nor the relevant dimensions, we may be hopeful 
that it does not lead to a different universality class. 

Random fields and random bonds: Uncorrelated fluctuations in the nearest neighbor 
coupling strengths (random bonds) in the presence of random field disorder do not break 
any new symmetries. Our random field Ising model fulfills two Harris criteria u/j35 > 2/d 
and v > 2/d. Adding random bond disorder cannot destroy the fixed point in the Harris 
criterion sense through added statistical fluctuations, because the random field disorder has 
already broken the relevant (translational) symmetry. It then seems plausible that systems 
with random bonds and random fields are in the same universality class as systems with 
random fields only. The ultimate justification for this conjecture may be drawn from the 
renormalization group picture. If the change in the generating functional due to the added 
new disorder turns out to be irrelevant under coarse graining, it will not affect the critical 
behavior on long length scales. Some preliminary studies seem to indicate that this would 
be case for random bonds in the presence of random fields. Further elaborations on this 
issue will be presented elsewhere [128]. 

Random bonds only: Similarly one might expect systems with random bonds only to 
be in the same universality class also. Because the critical magnetization M c = M(H C (R C )) 
is nonzero, the time reversal invariance will be broken at the critical point, just as it is 
broken in the case of random fields. The symmetries of the random field model and the 
random bond model would then be the same. In a soft spin model the dynamics could also 
be defined in the same way, using relaxational dynamics. One would then expect to see the 
same critical behavior on long length scales. In fact, in the random bond model one may 
consider the spins that flip outside the critical region to act as random fields for the spins 
that participate in the large avalanches near the critical point. We have already suggested 
that random bonds in the presence of random fields do not change the critical behavior. It 
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then seems plausible that the random bond problem would be in the same universality class 
also, as numerical simulations seem to confirm [50]. As a warning to the enthusiastic reader 
we should mention that the simplest mean field theory for this problem has some sicknesses 
that make it difficult to verify this conjecture in the RG framework. We believe, however, 
that these subtleties occur only in the infinite range mean-field theory, and that they are 
unimportant for the behavior in finite dimensions. Some simulations of the infinite range 
model with random bonds are reported in [118,116]. 

Random anisotropies: Realistic models of Barkhausen noise in polycrystalline mag- 
nets usually involve random anisotropies rather than random fields. In the same way it 
appears plausible that a nonequilibrium 0{n) model with random anisotropies [7,35,123] 
may be in the same universality class as the nonequilibrium RFIM. The external magnetic 
driving field breaks the rotational symmetry and time reversal invariance, just as in the case 
of random bonds. Again, spins that do not flip in the critical region may act as random fields 
for the spins participating in avalanches near the critical point, so that the essential features 
are the same as in our model, and one may expect to see the same critical exponents. 28 
The 0{n) model with random anisotropies is very similar to a continuous scalar spin model 
with random couplings to the external magnetic field (random "g- factors" ) . The mean- field 
theory for the random g-factor model turns out to have the same critical exponents as our 
random-field Ising model. There are no new terms generated in the RG description of this 
model either, it is therefore expected to be in the same universality class as our model. By 
symmetry we would neither expect any change in the exponents if there was randomness 
added through a distribution in the soft-spin potential well curvatures k (see our definition 



28 It may be that in some strong coupling limit the system will lose the ability to avalanche and 
all spins will smoothly rotate from down to up as the external magnetic field is increased. Our 
discussion above refers to the case where the coupling is weaker and avalanches do occur, as of 
course they do experimentally. 
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of the soft-spin potential V(sj) in eq. (31)), nor if random bonds are added to the system, 
as may be the case in real experimental systems. 

The RG formalism developed in this paper can be used as a convenient tool to verify 
these conjectures. One can write down the most general generating functional and verify 
for each of these models whether on long length scales the same kinds of terms become 
important or irrelevant as in our model. 

Long range interactions: 

The question about the effect of long range interactions is of equal importance. Depend- 
ing on the sample shape, dipole-dipole interactions can lead to long-range, antiferromagnetic 
interaction forces which are the reason for the breakup of the magnetization into Weiss- 
domains in conventional magnets [35,7]. In the case of martensites there are long range 
antiferroelastic strain fields present [58,8]. In references [126,129] we give an example of a 
critical exponent in a system with long range forces (from avalanche duration measurements 
in martensites [13]) that appears to be quite different from the corresponding exponent in 
our model. On the other hand, measurements of Barkhausen-noise distributions in magnets 
in the presence of long range demagnetizing fields seem to yield a critical exponent quite 
close to the corresponding exponent in our model [97]. 

In a recent preprint [97] Urbach, Madison and Markert study a model for a single moving 
domain wall without overhangs in the presence of infinite range antiferromagnetic interac- 
tions and quenched (random field) disorder. In an infinite system their model self-organizes 29 
without necessary parameter tuning to the same critical state seen in the absence of the in- 
finite range interactions right at the interface depinning threshold [28]. An analysis [124] of 
our ferromagnetic RFIM in the presence of infinite range antiferromagnetic interactions leads 
to an unchanged critical behavior except for a tilt of the entire magnetization curve in the 



This self-organization to the critical point is similar to the trivial self-organization expected in 
an experiment in the presence of a gradient field, which we proposed in references [126,129] 
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(M, H) plane: here too it does not change the critical properties. It would be interesting to 
see how these results would change for more physical long-range interactions. Dipole-dipole 
interactions decaying with distance as 1/x 3 , for example, might be appropriate. 

3. Thermal fluctuations 
(a) The equilibrium random field Ising model: 

The equilibrium properties of the random field Ising model, in particular the phase 
transition from paramagnetic to ferromagnetic (long range ordered) behavior, have been the 
subject of much controversy since the 1970s [62]. The reason is intriguing: experimental and 
theoretical studies of the approach to equilibrium show that near the critical temperature 
there seems to appear a "glassy" regime where relaxation to equilibrium becomes very slow. 
Activated by thermal fluctuations the system tumbles over free energy barriers to lower and 
lower valleys in the free energy landscape, until it has reached the lowest possible state, the 
equilibrium or ground state. The higher those barriers are compared to the typical energy 
of thermal fluctuations, the longer the relaxation process takes. At low temperatures, due 
to the effect of disorder, some of these barriers are so large (diverging in an infinite system), 
that the system gets stuck in some metastable state and never reaches true equilibrium on 
measurement time scales. On long length scales (and experimental time scales) thermal 
fluctuations become irrelevant and collective behavior emerges. When driven by an external 
field, the system moves through a local valley in the free energy landscape, and collective 
behavior in the form of avalanches is found when the system reaches a descending slope 
in the free energy surface. The present state of the system depends on its history — a 
phenomenon commonly observed as hysteresis, 
(b) The nonequilibrium random field Ising model: 

We have studied this hysteresis in the zero temperature random field Ising model, far 
from equilibrium and in the absence of any thermal fluctuations. We found a critical point, 
at which the shape of the hysteresis loops (magnetization versus magnetic field) changes 
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continuously from displaying a jump in the magnetization to a smooth curve. Interestingly 
the nonequilibrium critical exponents associated with the universal behavior near this point 
in d = 3 dimensions seem to match those obtained from 3 dimensional simulations of the 
equilibrium phase transition point approximately within the error bars (see table I). This is 
surprising, since the physical starting points of the two systems are very different. Further- 
more, our perturbation expansion in e = 6 — d for nonequilibrium critical exponents can be 
mapped onto the expansion for the equilibrium problem to all orders in e. Our expansion 
stems from a dynamical systems description of a deterministic process, which takes into 
account the history of the system and is designed to single out the correct metastable state, 
while the calculation for the equilibrium problem involves temperature fluctuations and no 
history dependence at all. 

(c) The crossover: 

It would be interesting to see if there is actually a deeper connection between the nonequi- 
librium and equilibrium critical points, and whether the calculation for the nonequilibrium 
model could be used to resolve long-standing difficulties with the perturbation expansion for 
the equilibrium model. The idea is to introduce temperature fluctuations in the nonequilib- 
rium calculation, and at the same time a finite sweeping frequency for the external driving 
force. The lower the sweeping frequency Q at fixed temperature, the more equilibrated 
the system and the longer the length scale above which nonequilibrium behavior emerges. 
Tuning fl would allow one to explore the whole crossover region between the two extreme 
cases that are found in the literature (far from and close to equilibrium). Contrary to pre- 
vious treatments of relaxation, the history dependence that is so essential in experimental 
realizations, emerges naturally from this approach. At fixed temperature, but for progres- 
sively lower sweeping frequencies, one expects to see smaller and smaller hysteresis loops, 
asymptotically attaining a universal shape at low enough frequencies. The tails of these 
hysteresis loops will match the equilibrium magnetization curve. In the limit of zero fre- 
quency, the hysteresis loop shrinks to a point, and equilibrium is expected at all values of the 



116 



external magnetic field. On the other hand, taking temperature to zero first, should yield 
nonequilibrium behavior as seen in our recent work. The prospect of relating equilibrium 
and nonequilibrium critical behavior as two limits at opposite edges of the experimentally 
relevant crossover regime is an exciting challenge. 
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FIGURES 



FIG. 1. Experiment: Magnetic hysteresis loops of a 60 nm thick Gd film for various 
annealing temperatures (as indicated next to each loop) and constant annealing time (3 minutes 
each). All measurements are performed at 200 ± 5K. The sweeping frequency of the external 
magnetic field is 0.5 Hz (from A. Berger, unpublished). 

FIG. 2. Experiment: Distribution of pulse areas (p), integrated over the hysteresis 
loop for 81% Ni-Fe wires after various heat treatments. The originally hard drawn wires 
have been subjected to a one- hour heat treatment in high vacuum at temperatures of 240° C or 
460° (7 and cooled down in the furnace. (From U. Lieneweg and W. Grosse-Nobis, Int. J. Magn 3, 
11 (1972).) 

FIG. 3. Equilibrium magnetization curve M(H) for the pure Ising model at zero temper- 
ature. 

FIG. 4. Nonequilibrium magnetization curve M(H) in the pure Ising model at zero tem- 
perature for the dynamics defined in the text. 

FIG. 5. Mean-field magnetization curves for the nonequilibrium zero temperature random 
field Ising model at various values of the disorder R = 0.6 J < R c (a), R = R C = ^Jpi^J = 0.798 J 
(b), and R = J > R c (c). 

FIG. 6. Mean-field magnetization curves for the soft-spin version of the zero tem- 
perature random field Ising model at various values of the disorder R = 1.3J < R c (a), 
R = R C = 2kJ/({k - J)v^7r) = 1.6 J (see appendix A 6) (b), and R = 2 J > R c (c). 

FIG. 7. Mean-field phase diagram for the nonequilibrium zero temperature random field 
Ising model. The critical point studied in this paper is at R = R c , H = H C (R C ), with H C (R C ) = 
in the hard-spin mean field theory. There are two relevant directions r = {R c — R)/R and 
h = H — H C (R C ) near this critical point. The bold line indicates the threshold field H™(R) for 
the onset of the infinite avalanche upon monotonically increasing the external magnetic field. The 
dashed line describes H l c (R) for a decreasing external magnetic field. The three dotted vertical lines 
marked (a), (b), and (c) describe the paths in parameter space which lead to the corresponding 
hysteresis loops shown in figure 5. 
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FIG. 8. Mean-field phase diagram for the soft-spin version of the nonequilibrium zero 
temperature random field Ising model. The diagram is plotted analogously to figure 7. Magnetic 
field sweeps along the lines (a), (b) and (c) lead to the corresponding soft-spin hysteresis curves 
shown in figure 6. Note that here, in contrast to the hard-spin model the value of the critical 
field, H C (R C ) does depend on the history of the system: for monotonically increasing external 
magnetic field H C (R C ) = H^(R C ) = k — J, and for monotonically decreasing external magnetic field 
H C {R C ) = H l c (R c ) = —(k — J) (see appendix A 6). This implies that in contrast to the hard-spin 
mean-field theory of figure 7, the soft-spin mean-field theory displays hysteresis for all finite disorder 
values, i.e. even at R > R c . 

FIG. 9. Simulated hysteresis curves for two small realizations of the nonequilibrium 
RFIM in three dimensions. Each sample consists of only 5 3 spins, with periodic boundary 
conditions. The two systems have different configurations of random fields that are taken from the 
same distribution p(f) with standard deviation R = 5 J > R c . (In 3 dimensions R c = (2. 16 ±0.03) J 
[45].) Note that here, as in all plots of numerical simulation results in finite dimensions J denotes 
the strength of the nearest neighbor coupling J^-, which differs from its definition in the analytical 
calculation by the coordination coordination number z - see footnote 4 on page 16. 

FIG. 10. Mean-field avalanche size distribution integrated over the hysteresis loop for 
systems with 1000000 spins at various disorder values R > R c = 0.798 J: (a) R = 1.46J (averaged 
over 10 different configurations of random fields), (b) R = 1.069J (averaged over 5 different config- 
urations of random fields), and (c) R = 0.912 J (averaged over 10 different configurations of random 
fields). Each curve is a histogram of all avalanche sizes found as the magnetic field is raised from 
— oo to +oo, normalized by the number of spins in the system. For small \r\ = \R C — R\/R the distri- 
bution roughly follows a power law D(S,r) ~ ( 5 , -( t +°"p < ') U p to a certain cutoff size Syjidx ^1^1 
which scales to infinity as r is taken to zero. The straight line above the three data curves in the 
figure represents an extrapolation to the critical point R = R c in an infinite system, where one 
expects to see a pure power law distribution on all length scales D(S,r) ~ g-ir+afiS) wfth ^ e 
mean field values of the corresponding exponents r + a(35 = 2.25. 

FIG. 11. Feynman diagrams. The relevant corrections to first order in e = 6 — d for the 
constant part \~ l I J m the propagator (a), and for the vertex u (b). Figure (c) shows an example 
of a diagram forbidden by causality. 
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FIG. 12. Phase diagram and flows (schematic), (a) The vertical axis is the external 
field H, responsible for pulling the system from down to up. The horizontal axis is the width of 
the random-field distribution R. The bold line is H C (R), the location of the infinite avalanche 
(assuming an initial condition with all spins down and a slowly increasing external field). The 
critical point we study is the end point of the infinite avalanche line (R C ,H C (R C )). 

Using the analogy with the Ising model (see text) we also show the RG flows around the 
critical point. Here we ignore the RG motion of the critical point itself: equivalently, the figure 
can represent a section through the critical fixed point tangent to the two unstable eigenvectors 
(labeled h and r). Two systems on the same RG trajectory (dashed thin lines) have the same 
long-wavelength properties (correlation functions ...) except for an overall change in length scale, 
leading to the scaling collapse of equation (9). The r eigendirection to the left extends along the 
infinite avalanche line; to the right, we speculate that it lies along the percolation threshold for up 
spins. 

(b) 0(e) RG flows below 6 dimensions in the (x -1 ,^) plane (see text). Linearization around 
the Wilson-Fisher (WF) fixed point yields the exponents given to 0(e) in the table. In the vicinity 
of the repulsive u = = x^ 1 (MFT) fixed points one obtains the old mean- field exponents. 

FIG. 13. Feynman diagram. The correction to 0(e) to the vertex w in an expansion about 
8 dimensions, see eq. (110) in the text. 

FIG. 14. Borel resummed critical exponents and simulation results. Shown are the 
numerical values of the exponents 1/v, r/, and (35 = v(d — rj)/2 (triangles, diamonds, and circles 
respectively) in 3, 4, and 5 dimensions and in mean field theory (dimension 6 and higher). The 
error bars denote systematic errors in finding the exponents from collapses of curves at different 
values of disorder R. Statistical errors are smaller. The dashed lines are the Borel sums to fifth 
order in e for the same exponents (see text). 

FIG. 15. Comparison to numerical results. Numerical values (filled symbols) of the expo- 
nents T + a(35, t, 1/V, avz, and av (circles, diamond, triangles up, squares, and triangle left) in 2, 
3,4, and 5 dimensions. The empty symbols are values for these exponents in mean field (dimension 
6). Note that the value of r in 2d was not measured. The empty diamond represents the expected 
value [98,129]. The numerical results are courtesy of Olga Perkovic [45,98] from simulations of 
sizes up to 7000 2 , 1000 3 , 80 4 , and 50 5 spins, where for 320 3 for example, more than 700 different 
random field configurations were measured. The long-dashed lines are the e expansions to first 
order for the exponents r + a(38, r, avz, and av. They are: r + a(35 = f — f, t = § + 0(e 2 ), and 
avz = 7} + 0(e 2 ), and av = \ + 0(e 2 ) where e = 6 — d and d is the dimension. The short-dashed line 
is the Borel sum for 1/v to fifth order in e. The other exponents can be obtained from exponent 
equalities (see section IV I in the text). The error bars denote systematic errors in finding the 
exponents from collapses of curves at different values of disorder R. Statistical errors are smaller. 

FIG. 16. Contour lines for the correlation length in the (r',h r ) plane (schematic). The 
tilted coordinate axes indicate the physical directions (r, h). Since v/([35) < v, the correlation 
length £ changes faster in the (0, h') direction than in the (r',0) direction. 
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FIG. 17. Feynman diagrams. The perturbative expansion about mean-field theory is pre- 
sented here by Feynman diagrams, (a) Graph for the vertex u. Incoming arrows denote ij fields, 
outgoing arrows denote rj fields, (b) Example of a diagram which violates causality and is therefore 
forbidden, (c) Graph for the vertex U2,o- (d) Example of a diagram that is zero due to momentum 
conservation [78]. 

FIG. 18. Feynman diagram for the lowest order correction to U2fi(H\, H2). The magnetic 
fields corresponding to the times at which the vertices are evaluated are indicated. The propagators 
do not couple different fields. 

FIG. 19. Different Borel-resummations for v. Borel resummation of the perturbation 
series for v to 0(e 5 ) in 6 — e dimensions using (a) the method of Vladimirov, Kazakov, and Tarasov 
[106], which does not impose a pole for v or any other independent information, and (b) the results 
of Le Guillou and Zinn- Justin [110], which are obtained by explicitly assuming a singularity for 
v at a (variationally determined) critical value of e and by imposing the exactly known value in 
two dimensions. Curve (b) is based on an old result for the epsilon-expansion which later turned 
out to out to have the wrong 5th order term [59]. Our own Borel-resummation (curve (a)) on 
the other hand has been obtained using the newest, presumably correct result for the fifth order 
term [59]. Partly by design, curve (b) agrees very well with the value of v in the equilibrium pure 
Ising model in d — 2 dimensions. Le Guillou and Zinn- Justin quote an error due to truncation of 
the series, which is increasingly larger than 10% below 2.5 dimensions. Curve (a) agrees better 
with the numerical results in our zero-temperature avalanche model, indicated for the respective 
dimensions by the black diamonds with error bars [45]. It is believed that the two resummation 
methods should converge to the same results if taken to high enough order, though this has never 
been proven. Also shown are values for v for the equilibrium random field Ising model in three 
dimensions (circles) from different sources [87,88]. 

FIG. 20. Feynman diagram. The correction to 0(e) to the vertex w in an expansion about 
8 dimensions, see eq. (E3) in the text. 

FIG. 21. The function defined in equation (F10), plotted as a function of H±. In the 

figure, dH denotes the amplitude which is called AH in the text. 
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TABLES 



TABLE I. Numerical results for the critical exponents in three dimensions for our hys- 
teresis model [45,98] and for the equilibrium zero temperature random field Ising model [62,87,88]. 
The breakdown of hyperscaling exponent 9 is calculated for the the hysteresis model from the 
relation (3 + (35 = (d — 9)v (see section a and references [126,129]). The values of the critical 
exponents of the two models remain within each other's errorbars (except for v and perhaps r], 
although Dayan, Schwartz, and Young [89] found that v ~ 1.4 in the three dimensional equilibrium 
random-field Ising model from real space renormalization group calculations); the equality of the 
exponents was conjectured by Maritan et al [86]. The numerical agreement may not be so surpris- 
ing, if one remembers that the 6 — e expansion is the same for all exponents of the two models. 
Nevertheless there is always room for nonperturbative corrections, so that the exponents might 
still be different in 3 dimensions (see section VIII B, and figure 19). Physically, the agreement is 
rather unexpected, since the nature of the two models is very different. While the hysteresis model 
is far from equilibrium, occupying a history dependent, metastable state, the equilibrium RFIM is 
always in the lowest free energy state. One may speculate, however, about a presumably universal 
crossover from our hysteresis model to the equilibrium random field Ising model as temperature 
fluctuations and a finite field-sweeping frequency Q are introduced (see appendix G3). 





Hysteresis loop [45] 


Equilibrium RFIM [87,88] 


exponents 


in 3 dimensions 


in 3 dimensions 




(courtesy Olga Perkovic) 




V 


1.42 ±0.17 


0.97, 1.30, 1.02 ±0.06 


P 


0.0 ±0.43 


-0.1, 0.05, 0.06 ±0.07 


(55 


1.81 ±0.36 


1.6, 1.9 ±0.4, 1.83 ±0.18 


7] 


0.79 ±0.29 


0.25, 0.5 ±0.5, 0.14 ±0.067 


6 


1.5±0.5 


1.45, 1.5 ±0.45, 1.851 ±0.067 


R c (Gaussian) 


2.16 ±0.03 


2.3 ± 0.2 [88] 


H C (R C ) 


1.435 ±0.004 


(by symmetry) 
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Phase Diagram in Mean-Field Theory 
(hard- spin model) 
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Phase Diagram in Mean-Field Theory 
(soft- spin model) 

H 



A 



k, 






»k-J 
















< 1 


^~ 








(a) 


(b) 


(c) 



R 




10.0 



Magnetic Field H/J 




Avalanche size (S) 
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